| 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 // 26 //
27 // ------------------------------------------- 27 // --------------------------------------------------------------
28 // 28 //
29 // File name: G4LowEnergyIonisation 29 // File name: G4LowEnergyIonisation
30 // 30 //
31 // Author: Alessandra Forti, Vladimir I 31 // Author: Alessandra Forti, Vladimir Ivanchenko
32 // 32 //
33 // Creation date: March 1999 33 // Creation date: March 1999
34 // 34 //
35 // Modifications: 35 // Modifications:
36 // - 11.04.2000 VL 36 // - 11.04.2000 VL
37 // Changing use of float and G4float casts t 37 // Changing use of float and G4float casts to G4double casts
38 // because of problems with optimisation (bu 38 // because of problems with optimisation (bug ?)
39 // 10.04.2000 VL 39 // 10.04.2000 VL
40 // - Correcting Fluorescence transition probab 40 // - Correcting Fluorescence transition probabilities in order to take into account
41 // non-radiative transitions. No Auger elect 41 // non-radiative transitions. No Auger electron simulated yet: energy is locally deposited.
42 // 10.04.2000 VL 42 // 10.04.2000 VL
43 // - Correction of incident electron final mom 43 // - Correction of incident electron final momentum direction
44 // 07.04.2000 VL+LU 44 // 07.04.2000 VL+LU
45 // - First implementation of continuous energy 45 // - First implementation of continuous energy loss
46 // 22.03.2000 VL 46 // 22.03.2000 VL
47 // - 1 bug corrected in SelectRandomAtom metho 47 // - 1 bug corrected in SelectRandomAtom method (units)
48 // 17.02.2000 Veronique Lefebure 48 // 17.02.2000 Veronique Lefebure
49 // - 5 bugs corrected: 49 // - 5 bugs corrected:
50 // *in Fluorescence, 2 bugs affecting 50 // *in Fluorescence, 2 bugs affecting
51 // . localEnergyDeposition and 51 // . localEnergyDeposition and
52 // . number of emitted photons that was then 52 // . number of emitted photons that was then always 1 less
53 // *in EnergySampling method: 53 // *in EnergySampling method:
54 // . expon = Parms[13]+1; (instead of uncorr 54 // . expon = Parms[13]+1; (instead of uncorrect -1)
55 // . rejection /= Parms[6];(instead of uncor 55 // . rejection /= Parms[6];(instead of uncorrect Parms[7])
56 // . Parms[6] is apparently corrupted in the 56 // . Parms[6] is apparently corrupted in the data file (often = 0)
57 // -->Compute normalisation into local var 57 // -->Compute normalisation into local variable rejectionMax
58 // and use rejectionMax in stead of Parms 58 // and use rejectionMax in stead of Parms[6]
59 // 59 //
60 // Added Livermore data table construction met 60 // Added Livermore data table construction methods A. Forti
61 // Modified BuildMeanFreePath to read new data 61 // Modified BuildMeanFreePath to read new data tables A. Forti
62 // Added EnergySampling method A. Forti 62 // Added EnergySampling method A. Forti
63 // Modified PostStepDoIt to insert sampling wi 63 // Modified PostStepDoIt to insert sampling with EEDL data A. Forti
64 // Added SelectRandomAtom A. Forti 64 // Added SelectRandomAtom A. Forti
65 // Added map of the elements A. Forti 65 // Added map of the elements A. Forti
66 // 20.09.00 V.Ivanchenko update fluctuations 66 // 20.09.00 V.Ivanchenko update fluctuations
67 // 24.04.01 V.Ivanchenko remove RogueWave 67 // 24.04.01 V.Ivanchenko remove RogueWave
68 // 22.05.01 V.Ivanchenko update calculation of 68 // 22.05.01 V.Ivanchenko update calculation of delta-ray kinematic +
69 // clean up the code 69 // clean up the code
70 // 02.08.01 V.Ivanchenko fix energy conservati 70 // 02.08.01 V.Ivanchenko fix energy conservation for small steps
71 // 18.08.01 V.Ivanchenko fix energy conservati 71 // 18.08.01 V.Ivanchenko fix energy conservation for pathalogical delta-energy
72 // 01.10.01 E. Guardincerri Replaced fluoresce 72 // 01.10.01 E. Guardincerri Replaced fluorescence generation in PostStepDoIt
73 // according to desig 73 // according to design iteration
74 // 04.10.01 MGP Minor clean-up in 74 // 04.10.01 MGP Minor clean-up in the fluo section, removal of
75 // compilation warnin 75 // compilation warnings and extra protection to
76 // prevent from acces 76 // prevent from accessing a null pointer
77 // 29.09.01 V.Ivanchenko revision based on 77 // 29.09.01 V.Ivanchenko revision based on design iteration
78 // 10.10.01 MGP Revision to improv 78 // 10.10.01 MGP Revision to improve code quality and
79 // consistency with d 79 // consistency with design
80 // 18.10.01 V.Ivanchenko Add fluorescence A 80 // 18.10.01 V.Ivanchenko Add fluorescence AlongStepDoIt
81 // 18.10.01 MGP Revision to improv 81 // 18.10.01 MGP Revision to improve code quality and
82 // consistency with d 82 // consistency with design
83 // 19.10.01 V.Ivanchenko update according t 83 // 19.10.01 V.Ivanchenko update according to new design, V.Ivanchenko
84 // 26.10.01 V.Ivanchenko clean up deexcitat 84 // 26.10.01 V.Ivanchenko clean up deexcitation
85 // 28.10.01 V.Ivanchenko update printout 85 // 28.10.01 V.Ivanchenko update printout
86 // 29.11.01 V.Ivanchenko New parametrisatio 86 // 29.11.01 V.Ivanchenko New parametrisation introduced
87 // 25.03.02 V.Ivanchneko Fix in fluorescenc 87 // 25.03.02 V.Ivanchneko Fix in fluorescence
88 // 28.03.02 V.Ivanchenko Add flag of fluore 88 // 28.03.02 V.Ivanchenko Add flag of fluorescence
89 // 28.05.02 V.Ivanchenko Remove flag fStopA 89 // 28.05.02 V.Ivanchenko Remove flag fStopAndKill
90 // 31.05.02 V.Ivanchenko Add path of Fluo + 90 // 31.05.02 V.Ivanchenko Add path of Fluo + Auger cuts to
91 // AtomicDeexcitation 91 // AtomicDeexcitation
92 // 03.06.02 MGP Restore fStopAndKi 92 // 03.06.02 MGP Restore fStopAndKill
93 // 19.06.02 VI Additional printou 93 // 19.06.02 VI Additional printout
94 // 30.07.02 VI Fix in restricted 94 // 30.07.02 VI Fix in restricted energy loss
95 // 20.09.02 VI Remove ActivateFlu 95 // 20.09.02 VI Remove ActivateFlurescence from SetCut...
96 // 21.01.03 VI Cut per region 96 // 21.01.03 VI Cut per region
97 // 12.02.03 VI Change signature f 97 // 12.02.03 VI Change signature for Deexcitation
98 // 12.04.03 V.Ivanchenko Cut per region for 98 // 12.04.03 V.Ivanchenko Cut per region for fluo AlongStep
99 // 31.08.04 V.Ivanchenko Add density correc 99 // 31.08.04 V.Ivanchenko Add density correction
100 // 100 //
101 // ------------------------------------------- 101 // --------------------------------------------------------------
102 102
103 #include "G4LowEnergyIonisation.hh" 103 #include "G4LowEnergyIonisation.hh"
104 #include "G4PhysicalConstants.hh" 104 #include "G4PhysicalConstants.hh"
105 #include "G4SystemOfUnits.hh" 105 #include "G4SystemOfUnits.hh"
106 #include "G4RDeIonisationSpectrum.hh" 106 #include "G4RDeIonisationSpectrum.hh"
107 #include "G4RDeIonisationCrossSectionHandler.h 107 #include "G4RDeIonisationCrossSectionHandler.hh"
108 #include "G4RDAtomicTransitionManager.hh" 108 #include "G4RDAtomicTransitionManager.hh"
109 #include "G4RDAtomicShell.hh" 109 #include "G4RDAtomicShell.hh"
110 #include "G4RDVDataSetAlgorithm.hh" 110 #include "G4RDVDataSetAlgorithm.hh"
111 #include "G4RDSemiLogInterpolation.hh" 111 #include "G4RDSemiLogInterpolation.hh"
112 #include "G4RDLogLogInterpolation.hh" 112 #include "G4RDLogLogInterpolation.hh"
113 #include "G4RDEMDataSet.hh" 113 #include "G4RDEMDataSet.hh"
114 #include "G4RDVEMDataSet.hh" 114 #include "G4RDVEMDataSet.hh"
115 #include "G4RDCompositeEMDataSet.hh" 115 #include "G4RDCompositeEMDataSet.hh"
116 #include "G4EnergyLossTables.hh" 116 #include "G4EnergyLossTables.hh"
117 #include "G4RDShellVacancy.hh" 117 #include "G4RDShellVacancy.hh"
118 #include "G4UnitsTable.hh" 118 #include "G4UnitsTable.hh"
119 #include "G4Electron.hh" 119 #include "G4Electron.hh"
120 #include "G4Gamma.hh" 120 #include "G4Gamma.hh"
121 #include "G4ProductionCutsTable.hh" 121 #include "G4ProductionCutsTable.hh"
122 122
123 G4LowEnergyIonisation::G4LowEnergyIonisation(c 123 G4LowEnergyIonisation::G4LowEnergyIonisation(const G4String& nam)
124 : G4eLowEnergyLoss(nam), 124 : G4eLowEnergyLoss(nam),
125 crossSectionHandler(0), 125 crossSectionHandler(0),
126 theMeanFreePath(0), 126 theMeanFreePath(0),
127 energySpectrum(0), 127 energySpectrum(0),
128 shellVacancy(0) 128 shellVacancy(0)
129 { 129 {
130 cutForPhotons = 250.0*eV; 130 cutForPhotons = 250.0*eV;
131 cutForElectrons = 250.0*eV; 131 cutForElectrons = 250.0*eV;
132 verboseLevel = 0; 132 verboseLevel = 0;
133 } 133 }
134 134
135 135
136 G4LowEnergyIonisation::~G4LowEnergyIonisation( 136 G4LowEnergyIonisation::~G4LowEnergyIonisation()
137 { 137 {
138 delete crossSectionHandler; 138 delete crossSectionHandler;
139 delete energySpectrum; 139 delete energySpectrum;
140 delete theMeanFreePath; 140 delete theMeanFreePath;
141 delete shellVacancy; 141 delete shellVacancy;
142 } 142 }
143 143
144 144
145 void G4LowEnergyIonisation::BuildPhysicsTable( 145 void G4LowEnergyIonisation::BuildPhysicsTable(const G4ParticleDefinition& aParticleType)
146 { 146 {
147 if(verboseLevel > 0) { 147 if(verboseLevel > 0) {
148 G4cout << "G4LowEnergyIonisation::BuildPhy 148 G4cout << "G4LowEnergyIonisation::BuildPhysicsTable start"
149 << G4endl; 149 << G4endl;
150 } 150 }
151 151
152 cutForDelta.clear(); 152 cutForDelta.clear();
153 153
154 // Create and fill IonisationParameters once 154 // Create and fill IonisationParameters once
155 if( energySpectrum != 0 ) delete energySpect 155 if( energySpectrum != 0 ) delete energySpectrum;
156 energySpectrum = new G4RDeIonisationSpectrum 156 energySpectrum = new G4RDeIonisationSpectrum();
157 157
158 if(verboseLevel > 0) { 158 if(verboseLevel > 0) {
159 G4cout << "G4RDVEnergySpectrum is initiali 159 G4cout << "G4RDVEnergySpectrum is initialized"
160 << G4endl; 160 << G4endl;
161 } 161 }
162 162
163 // Create and fill G4RDCrossSectionHandler o 163 // Create and fill G4RDCrossSectionHandler once
164 164
165 if ( crossSectionHandler != 0 ) delete cross 165 if ( crossSectionHandler != 0 ) delete crossSectionHandler;
166 G4RDVDataSetAlgorithm* interpolation = new G 166 G4RDVDataSetAlgorithm* interpolation = new G4RDSemiLogInterpolation();
167 G4double lowKineticEnergy = GetLowerBoundEl 167 G4double lowKineticEnergy = GetLowerBoundEloss();
168 G4double highKineticEnergy = GetUpperBoundEl 168 G4double highKineticEnergy = GetUpperBoundEloss();
169 G4int totBin = GetNbinEloss(); 169 G4int totBin = GetNbinEloss();
170 crossSectionHandler = new G4RDeIonisationCro 170 crossSectionHandler = new G4RDeIonisationCrossSectionHandler(energySpectrum,
171 interpolation, 171 interpolation,
172 lowKineticEnergy, 172 lowKineticEnergy,
173 highKineticEnergy, 173 highKineticEnergy,
174 totBin); 174 totBin);
175 crossSectionHandler->LoadShellData("ioni/ion 175 crossSectionHandler->LoadShellData("ioni/ion-ss-cs-");
176 176
177 if (verboseLevel > 0) { 177 if (verboseLevel > 0) {
178 G4cout << GetProcessName() 178 G4cout << GetProcessName()
179 << " is created; Cross section data 179 << " is created; Cross section data: "
180 << G4endl; 180 << G4endl;
181 crossSectionHandler->PrintData(); 181 crossSectionHandler->PrintData();
182 G4cout << "Parameters: " 182 G4cout << "Parameters: "
183 << G4endl; 183 << G4endl;
184 energySpectrum->PrintData(); 184 energySpectrum->PrintData();
185 } 185 }
186 186
187 // Build loss table for IonisationIV 187 // Build loss table for IonisationIV
188 188
189 BuildLossTable(aParticleType); 189 BuildLossTable(aParticleType);
190 190
191 if(verboseLevel > 0) { 191 if(verboseLevel > 0) {
192 G4cout << "The loss table is built" 192 G4cout << "The loss table is built"
193 << G4endl; 193 << G4endl;
194 } 194 }
195 195
196 if (&aParticleType==G4Electron::Electron()) 196 if (&aParticleType==G4Electron::Electron()) {
197 197
198 RecorderOfElectronProcess[CounterOfElectro 198 RecorderOfElectronProcess[CounterOfElectronProcess] = (*this).theLossTable;
199 CounterOfElectronProcess++; 199 CounterOfElectronProcess++;
200 PrintInfoDefinition(); 200 PrintInfoDefinition();
201 201
202 } else { 202 } else {
203 203
204 RecorderOfPositronProcess[CounterOfPositro 204 RecorderOfPositronProcess[CounterOfPositronProcess] = (*this).theLossTable;
205 CounterOfPositronProcess++; 205 CounterOfPositronProcess++;
206 } 206 }
207 207
208 // Build mean free path data using cut value 208 // Build mean free path data using cut values
209 209
210 if( theMeanFreePath ) delete theMeanFreePath 210 if( theMeanFreePath ) delete theMeanFreePath;
211 theMeanFreePath = crossSectionHandler-> 211 theMeanFreePath = crossSectionHandler->
212 BuildMeanFreePathForMateri 212 BuildMeanFreePathForMaterials(&cutForDelta);
213 213
214 if(verboseLevel > 0) { 214 if(verboseLevel > 0) {
215 G4cout << "The MeanFreePath table is built 215 G4cout << "The MeanFreePath table is built"
216 << G4endl; 216 << G4endl;
217 if(verboseLevel > 1) theMeanFreePath->Prin 217 if(verboseLevel > 1) theMeanFreePath->PrintData();
218 } 218 }
219 219
220 // Build common DEDX table for all ionisatio 220 // Build common DEDX table for all ionisation processes
221 221
222 BuildDEDXTable(aParticleType); 222 BuildDEDXTable(aParticleType);
223 223
224 if (verboseLevel > 0) { 224 if (verboseLevel > 0) {
225 G4cout << "G4LowEnergyIonisation::BuildPhy 225 G4cout << "G4LowEnergyIonisation::BuildPhysicsTable end"
226 << G4endl; 226 << G4endl;
227 } 227 }
228 } 228 }
229 229
230 230
231 void G4LowEnergyIonisation::BuildLossTable(con 231 void G4LowEnergyIonisation::BuildLossTable(const G4ParticleDefinition& )
232 { 232 {
233 // Build table for energy loss due to soft b 233 // Build table for energy loss due to soft brems
234 // the tables are built for *MATERIALS* binn 234 // the tables are built for *MATERIALS* binning is taken from LowEnergyLoss
235 235
236 G4double lowKineticEnergy = GetLowerBoundEl 236 G4double lowKineticEnergy = GetLowerBoundEloss();
237 G4double highKineticEnergy = GetUpperBoundEl 237 G4double highKineticEnergy = GetUpperBoundEloss();
238 size_t totBin = GetNbinEloss(); 238 size_t totBin = GetNbinEloss();
239 239
240 // create table 240 // create table
241 241
242 if (theLossTable) { 242 if (theLossTable) {
243 theLossTable->clearAndDestroy(); 243 theLossTable->clearAndDestroy();
244 delete theLossTable; 244 delete theLossTable;
245 } 245 }
246 const G4ProductionCutsTable* theCoupleTable= 246 const G4ProductionCutsTable* theCoupleTable=
247 G4ProductionCutsTable::GetProductionCu 247 G4ProductionCutsTable::GetProductionCutsTable();
248 size_t numOfCouples = theCoupleTable->GetTab 248 size_t numOfCouples = theCoupleTable->GetTableSize();
249 theLossTable = new G4PhysicsTable(numOfCoupl 249 theLossTable = new G4PhysicsTable(numOfCouples);
250 250
251 if (shellVacancy != 0) delete shellVacancy; 251 if (shellVacancy != 0) delete shellVacancy;
252 shellVacancy = new G4RDShellVacancy(); 252 shellVacancy = new G4RDShellVacancy();
253 G4DataVector* ksi = 0; 253 G4DataVector* ksi = 0;
254 G4DataVector* energy = 0; 254 G4DataVector* energy = 0;
255 size_t binForFluo = totBin/10; 255 size_t binForFluo = totBin/10;
256 256
257 G4PhysicsLogVector* bVector = new G4PhysicsL 257 G4PhysicsLogVector* bVector = new G4PhysicsLogVector(lowKineticEnergy,
258 highKineticEner 258 highKineticEnergy,
259 binForFluo); 259 binForFluo);
260 const G4RDAtomicTransitionManager* transitio 260 const G4RDAtomicTransitionManager* transitionManager = G4RDAtomicTransitionManager::Instance();
261 261
262 // Clean up the vector of cuts 262 // Clean up the vector of cuts
263 263
264 cutForDelta.clear(); 264 cutForDelta.clear();
265 265
266 // Loop for materials 266 // Loop for materials
267 267
268 for (size_t m=0; m<numOfCouples; m++) { 268 for (size_t m=0; m<numOfCouples; m++) {
269 269
270 // create physics vector and fill it 270 // create physics vector and fill it
271 G4PhysicsLogVector* aVector = new G4Physic 271 G4PhysicsLogVector* aVector = new G4PhysicsLogVector(lowKineticEnergy,
272 highKineticEnergy, 272 highKineticEnergy,
273 totBin); 273 totBin);
274 274
275 // get material parameters needed for the 275 // get material parameters needed for the energy loss calculation
276 const G4MaterialCutsCouple* couple = theCo 276 const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(m);
277 const G4Material* material= couple->GetMat 277 const G4Material* material= couple->GetMaterial();
278 278
279 // the cut cannot be below lowest limit 279 // the cut cannot be below lowest limit
280 G4double tCut = (*(theCoupleTable->GetEner 280 G4double tCut = (*(theCoupleTable->GetEnergyCutsVector(1)))[m];
281 if(tCut > highKineticEnergy) tCut = highKi 281 if(tCut > highKineticEnergy) tCut = highKineticEnergy;
282 cutForDelta.push_back(tCut); 282 cutForDelta.push_back(tCut);
283 const G4ElementVector* theElementVector = 283 const G4ElementVector* theElementVector = material->GetElementVector();
284 size_t NumberOfElements = material->GetNum 284 size_t NumberOfElements = material->GetNumberOfElements() ;
285 const G4double* theAtomicNumDensityVector 285 const G4double* theAtomicNumDensityVector =
286 material->GetAtomicNumDens 286 material->GetAtomicNumDensityVector();
287 if(verboseLevel > 0) { 287 if(verboseLevel > 0) {
288 G4cout << "Energy loss for material # " 288 G4cout << "Energy loss for material # " << m
289 << " tCut(keV)= " << tCut/keV 289 << " tCut(keV)= " << tCut/keV
290 << G4endl; 290 << G4endl;
291 } 291 }
292 292
293 // now comes the loop for the kinetic ener 293 // now comes the loop for the kinetic energy values
294 for (size_t i = 0; i<totBin; i++) { 294 for (size_t i = 0; i<totBin; i++) {
295 295
296 G4double lowEdgeEnergy = aVector->GetLow 296 G4double lowEdgeEnergy = aVector->GetLowEdgeEnergy(i);
297 G4double ionloss = 0.; 297 G4double ionloss = 0.;
298 298
299 // loop for elements in the material 299 // loop for elements in the material
300 for (size_t iel=0; iel<NumberOfElements; 300 for (size_t iel=0; iel<NumberOfElements; iel++ ) {
301 301
302 G4int Z = (G4int)((*theElementVector)[ 302 G4int Z = (G4int)((*theElementVector)[iel]->GetZ());
303 303
304 G4int nShells = transitionManager->NumberOfS 304 G4int nShells = transitionManager->NumberOfShells(Z);
305 305
306 for (G4int n=0; n<nShells; n++) { 306 for (G4int n=0; n<nShells; n++) {
307 307
308 G4double e = energySpectrum->Average 308 G4double e = energySpectrum->AverageEnergy(Z, 0.0, tCut,
309 309 lowEdgeEnergy, n);
310 G4double cs= crossSectionHandler->Fi 310 G4double cs= crossSectionHandler->FindValue(Z, lowEdgeEnergy, n);
311 ionloss += e * cs * theAtomicNumDe 311 ionloss += e * cs * theAtomicNumDensityVector[iel];
312 312
313 if(verboseLevel > 1 || (Z == 14 && l 313 if(verboseLevel > 1 || (Z == 14 && lowEdgeEnergy>1. && lowEdgeEnergy<0.)) {
314 G4cout << "Z= " << Z 314 G4cout << "Z= " << Z
315 << " shell= " << n 315 << " shell= " << n
316 << " E(keV)= " << lowEdgeEn 316 << " E(keV)= " << lowEdgeEnergy/keV
317 << " Eav(keV)= " << e/keV 317 << " Eav(keV)= " << e/keV
318 << " cs= " << cs 318 << " cs= " << cs
319 << " loss= " << ionloss 319 << " loss= " << ionloss
320 << " rho= " << theAtomicNumDensit 320 << " rho= " << theAtomicNumDensityVector[iel]
321 << G4endl; 321 << G4endl;
322 } 322 }
323 } 323 }
324 G4double esp = energySpectrum->Excitat 324 G4double esp = energySpectrum->Excitation(Z, lowEdgeEnergy);
325 ionloss += esp * theAtomicNumDensity 325 ionloss += esp * theAtomicNumDensityVector[iel];
326 326
327 } 327 }
328 if(verboseLevel > 1 || (m == 0 && lowEdg 328 if(verboseLevel > 1 || (m == 0 && lowEdgeEnergy>=1. && lowEdgeEnergy<=0.)) {
329 G4cout << "Sum: " 329 G4cout << "Sum: "
330 << " E(keV)= " << lowEdgeEn 330 << " E(keV)= " << lowEdgeEnergy/keV
331 << " loss(MeV/mm)= " << ionloss*m 331 << " loss(MeV/mm)= " << ionloss*mm/MeV
332 << G4endl; 332 << G4endl;
333 } 333 }
334 aVector->PutValue(i,ionloss); 334 aVector->PutValue(i,ionloss);
335 } 335 }
336 theLossTable->insert(aVector); 336 theLossTable->insert(aVector);
337 337
338 // fill data for fluorescence 338 // fill data for fluorescence
339 339
340 G4RDVDataSetAlgorithm* interp = new G4RDLo 340 G4RDVDataSetAlgorithm* interp = new G4RDLogLogInterpolation();
341 G4RDVEMDataSet* xsis = new G4RDCompositeEM 341 G4RDVEMDataSet* xsis = new G4RDCompositeEMDataSet(interp, 1., 1.);
342 for (size_t iel=0; iel<NumberOfElements; i 342 for (size_t iel=0; iel<NumberOfElements; iel++ ) {
343 343
344 G4int Z = (G4int)((*theElementVector)[ie 344 G4int Z = (G4int)((*theElementVector)[iel]->GetZ());
345 energy = new G4DataVector(); 345 energy = new G4DataVector();
346 ksi = new G4DataVector(); 346 ksi = new G4DataVector();
347 347
348 for (size_t j = 0; j<binForFluo; j++) { 348 for (size_t j = 0; j<binForFluo; j++) {
349 349
350 G4double lowEdgeEnergy = bVector->GetL 350 G4double lowEdgeEnergy = bVector->GetLowEdgeEnergy(j);
351 G4double cross = 0.; 351 G4double cross = 0.;
352 G4double eAverage= 0.; 352 G4double eAverage= 0.;
353 G4int nShells = transitionManager->NumberOfS 353 G4int nShells = transitionManager->NumberOfShells(Z);
354 354
355 for (G4int n=0; n<nShells; n++) { 355 for (G4int n=0; n<nShells; n++) {
356 356
357 G4double e = energySpectrum->Average 357 G4double e = energySpectrum->AverageEnergy(Z, 0.0, tCut,
358 358 lowEdgeEnergy, n);
359 G4double pro = energySpectrum->Proba 359 G4double pro = energySpectrum->Probability(Z, 0.0, tCut,
360 360 lowEdgeEnergy, n);
361 G4double cs= crossSectionHandler->Fi 361 G4double cs= crossSectionHandler->FindValue(Z, lowEdgeEnergy, n);
362 eAverage += e * cs * theAtomicNumD 362 eAverage += e * cs * theAtomicNumDensityVector[iel];
363 cross += cs * pro * theAtomicNu 363 cross += cs * pro * theAtomicNumDensityVector[iel];
364 if(verboseLevel > 1) { 364 if(verboseLevel > 1) {
365 G4cout << "Z= " << Z 365 G4cout << "Z= " << Z
366 << " shell= " << n 366 << " shell= " << n
367 << " E(keV)= " << lowEdgeEn 367 << " E(keV)= " << lowEdgeEnergy/keV
368 << " Eav(keV)= " << e/keV 368 << " Eav(keV)= " << e/keV
369 << " pro= " << pro 369 << " pro= " << pro
370 << " cs= " << cs 370 << " cs= " << cs
371 << G4endl; 371 << G4endl;
372 } 372 }
373 } 373 }
374 374
375 G4double coeff = 0.0; 375 G4double coeff = 0.0;
376 if(eAverage > 0.) { 376 if(eAverage > 0.) {
377 coeff = cross/eAverage; 377 coeff = cross/eAverage;
378 eAverage /= cross; 378 eAverage /= cross;
379 } 379 }
380 380
381 if(verboseLevel > 1) { 381 if(verboseLevel > 1) {
382 G4cout << "Ksi Coefficient for Z= 382 G4cout << "Ksi Coefficient for Z= " << Z
383 << " E(keV)= " << lowEdgeEn 383 << " E(keV)= " << lowEdgeEnergy/keV
384 << " Eav(keV)= " << eAverag 384 << " Eav(keV)= " << eAverage/keV
385 << " coeff= " << coeff 385 << " coeff= " << coeff
386 << G4endl; 386 << G4endl;
387 } 387 }
388 388
389 energy->push_back(lowEdgeEnergy); 389 energy->push_back(lowEdgeEnergy);
390 ksi->push_back(coeff); 390 ksi->push_back(coeff);
391 } 391 }
392 interp = new G4RDLogLogInterpolation(); 392 interp = new G4RDLogLogInterpolation();
393 G4RDVEMDataSet* set = new G4RDEMDataSet( 393 G4RDVEMDataSet* set = new G4RDEMDataSet(Z,energy,ksi,interp,1.,1.);
394 xsis->AddComponent(set); 394 xsis->AddComponent(set);
395 } 395 }
396 if(verboseLevel) xsis->PrintData(); 396 if(verboseLevel) xsis->PrintData();
397 shellVacancy->AddXsiTable(xsis); 397 shellVacancy->AddXsiTable(xsis);
398 } 398 }
399 delete bVector; 399 delete bVector;
400 } 400 }
401 401
402 402
403 G4VParticleChange* G4LowEnergyIonisation::Post 403 G4VParticleChange* G4LowEnergyIonisation::PostStepDoIt(const G4Track& track,
404 const G4Step& step) 404 const G4Step& step)
405 { 405 {
406 // Delta electron production mechanism on ba 406 // Delta electron production mechanism on base of the model
407 // J. Stepanek " A program to determine the 407 // J. Stepanek " A program to determine the radiation spectra due
408 // to a single atomic subshell ionisation by 408 // to a single atomic subshell ionisation by a particle or due to
409 // deexcitation or decay of radionuclides", 409 // deexcitation or decay of radionuclides",
410 // Comp. Phys. Comm. 1206 pp 1-19 (1997) 410 // Comp. Phys. Comm. 1206 pp 1-19 (1997)
411 411
412 aParticleChange.Initialize(track); 412 aParticleChange.Initialize(track);
413 413
414 const G4MaterialCutsCouple* couple = track.G 414 const G4MaterialCutsCouple* couple = track.GetMaterialCutsCouple();
415 G4double kineticEnergy = track.GetKineticEne 415 G4double kineticEnergy = track.GetKineticEnergy();
416 416
417 // Select atom and shell 417 // Select atom and shell
418 418
419 G4int Z = crossSectionHandler->SelectRandomA 419 G4int Z = crossSectionHandler->SelectRandomAtom(couple, kineticEnergy);
420 G4int shell = crossSectionHandler->SelectRan 420 G4int shell = crossSectionHandler->SelectRandomShell(Z, kineticEnergy);
421 const G4RDAtomicShell* atomicShell = 421 const G4RDAtomicShell* atomicShell =
422 (G4RDAtomicTransitionManager:: 422 (G4RDAtomicTransitionManager::Instance())->Shell(Z, shell);
423 G4double bindingEnergy = atomicShell->Bindin 423 G4double bindingEnergy = atomicShell->BindingEnergy();
424 G4int shellId = atomicShell->ShellId(); 424 G4int shellId = atomicShell->ShellId();
425 425
426 // Sample delta energy 426 // Sample delta energy
427 427
428 G4int index = couple->GetIndex(); 428 G4int index = couple->GetIndex();
429 G4double tCut = cutForDelta[index]; 429 G4double tCut = cutForDelta[index];
430 G4double tmax = energySpectrum->MaxEnergyO 430 G4double tmax = energySpectrum->MaxEnergyOfSecondaries(kineticEnergy);
431 G4double tDelta = energySpectrum->SampleEner 431 G4double tDelta = energySpectrum->SampleEnergy(Z, tCut, tmax,
432 432 kineticEnergy, shell);
433 433
434 if(tDelta == 0.0) 434 if(tDelta == 0.0)
435 return G4VContinuousDiscreteProcess::PostS 435 return G4VContinuousDiscreteProcess::PostStepDoIt(track, step);
436 436
437 // Transform to shell potential 437 // Transform to shell potential
438 G4double deltaKinE = tDelta + 2.0*bindingEne 438 G4double deltaKinE = tDelta + 2.0*bindingEnergy;
439 G4double primaryKinE = kineticEnergy + 2.0*b 439 G4double primaryKinE = kineticEnergy + 2.0*bindingEnergy;
440 440
441 // sampling of scattering angle neglecting a 441 // sampling of scattering angle neglecting atomic motion
442 G4double deltaMom = std::sqrt(deltaKinE*(del 442 G4double deltaMom = std::sqrt(deltaKinE*(deltaKinE + 2.0*electron_mass_c2));
443 G4double primaryMom = std::sqrt(primaryKinE* 443 G4double primaryMom = std::sqrt(primaryKinE*(primaryKinE + 2.0*electron_mass_c2));
444 444
445 G4double cost = deltaKinE * (primaryKinE + 2 445 G4double cost = deltaKinE * (primaryKinE + 2.0*electron_mass_c2)
446 / (deltaMom * prim 446 / (deltaMom * primaryMom);
447 447
448 if (cost > 1.) cost = 1.; 448 if (cost > 1.) cost = 1.;
449 G4double sint = std::sqrt(1. - cost*cost); 449 G4double sint = std::sqrt(1. - cost*cost);
450 G4double phi = twopi * G4UniformRand(); 450 G4double phi = twopi * G4UniformRand();
451 G4double dirx = sint * std::cos(phi); 451 G4double dirx = sint * std::cos(phi);
452 G4double diry = sint * std::sin(phi); 452 G4double diry = sint * std::sin(phi);
453 G4double dirz = cost; 453 G4double dirz = cost;
454 454
455 // Rotate to incident electron direction 455 // Rotate to incident electron direction
456 G4ThreeVector primaryDirection = track.GetMo 456 G4ThreeVector primaryDirection = track.GetMomentumDirection();
457 G4ThreeVector deltaDir(dirx,diry,dirz); 457 G4ThreeVector deltaDir(dirx,diry,dirz);
458 deltaDir.rotateUz(primaryDirection); 458 deltaDir.rotateUz(primaryDirection);
459 dirx = deltaDir.x(); 459 dirx = deltaDir.x();
460 diry = deltaDir.y(); 460 diry = deltaDir.y();
461 dirz = deltaDir.z(); 461 dirz = deltaDir.z();
462 462
463 463
464 // Take into account atomic motion del is re 464 // Take into account atomic motion del is relative momentum of the motion
465 // kinetic energy of the motion == bindingEn 465 // kinetic energy of the motion == bindingEnergy in V.Ivanchenko model
466 466
467 cost = 2.0*G4UniformRand() - 1.0; 467 cost = 2.0*G4UniformRand() - 1.0;
468 sint = std::sqrt(1. - cost*cost); 468 sint = std::sqrt(1. - cost*cost);
469 phi = twopi * G4UniformRand(); 469 phi = twopi * G4UniformRand();
470 G4double del = std::sqrt(bindingEnergy *(bin 470 G4double del = std::sqrt(bindingEnergy *(bindingEnergy + 2.0*electron_mass_c2))
471 / deltaMom; 471 / deltaMom;
472 dirx += del* sint * std::cos(phi); 472 dirx += del* sint * std::cos(phi);
473 diry += del* sint * std::sin(phi); 473 diry += del* sint * std::sin(phi);
474 dirz += del* cost; 474 dirz += del* cost;
475 475
476 // Find out new primary electron direction 476 // Find out new primary electron direction
477 G4double finalPx = primaryMom*primaryDirecti 477 G4double finalPx = primaryMom*primaryDirection.x() - deltaMom*dirx;
478 G4double finalPy = primaryMom*primaryDirecti 478 G4double finalPy = primaryMom*primaryDirection.y() - deltaMom*diry;
479 G4double finalPz = primaryMom*primaryDirecti 479 G4double finalPz = primaryMom*primaryDirection.z() - deltaMom*dirz;
480 480
481 // create G4DynamicParticle object for delta 481 // create G4DynamicParticle object for delta ray
482 G4DynamicParticle* theDeltaRay = new G4Dynam 482 G4DynamicParticle* theDeltaRay = new G4DynamicParticle();
483 theDeltaRay->SetKineticEnergy(tDelta); 483 theDeltaRay->SetKineticEnergy(tDelta);
484 G4double norm = 1.0/std::sqrt(dirx*dirx + di 484 G4double norm = 1.0/std::sqrt(dirx*dirx + diry*diry + dirz*dirz);
485 dirx *= norm; 485 dirx *= norm;
486 diry *= norm; 486 diry *= norm;
487 dirz *= norm; 487 dirz *= norm;
488 theDeltaRay->SetMomentumDirection(dirx, diry 488 theDeltaRay->SetMomentumDirection(dirx, diry, dirz);
489 theDeltaRay->SetDefinition(G4Electron::Elect 489 theDeltaRay->SetDefinition(G4Electron::Electron());
490 490
491 G4double theEnergyDeposit = bindingEnergy; 491 G4double theEnergyDeposit = bindingEnergy;
492 492
493 // fill ParticleChange 493 // fill ParticleChange
494 // changed energy and momentum of the actual 494 // changed energy and momentum of the actual particle
495 495
496 G4double finalKinEnergy = kineticEnergy - tD 496 G4double finalKinEnergy = kineticEnergy - tDelta - theEnergyDeposit;
497 if(finalKinEnergy < 0.0) { 497 if(finalKinEnergy < 0.0) {
498 theEnergyDeposit += finalKinEnergy; 498 theEnergyDeposit += finalKinEnergy;
499 finalKinEnergy = 0.0; 499 finalKinEnergy = 0.0;
500 aParticleChange.ProposeTrackStatus(fStopAn 500 aParticleChange.ProposeTrackStatus(fStopAndKill);
501 501
502 } else { 502 } else {
503 503
504 G4double norm = 1.0/std::sqrt(finalPx*fina 504 G4double norm = 1.0/std::sqrt(finalPx*finalPx+finalPy*finalPy+finalPz*finalPz);
505 finalPx *= norm; 505 finalPx *= norm;
506 finalPy *= norm; 506 finalPy *= norm;
507 finalPz *= norm; 507 finalPz *= norm;
508 aParticleChange.ProposeMomentumDirection(f 508 aParticleChange.ProposeMomentumDirection(finalPx, finalPy, finalPz);
509 } 509 }
510 510
511 aParticleChange.ProposeEnergy(finalKinEnergy 511 aParticleChange.ProposeEnergy(finalKinEnergy);
512 512
513 // Generation of Fluorescence and Auger 513 // Generation of Fluorescence and Auger
514 size_t nSecondaries = 0; 514 size_t nSecondaries = 0;
515 size_t totalNumber = 1; 515 size_t totalNumber = 1;
516 std::vector<G4DynamicParticle*>* secondaryVe 516 std::vector<G4DynamicParticle*>* secondaryVector = 0;
517 G4DynamicParticle* aSecondary = 0; 517 G4DynamicParticle* aSecondary = 0;
518 G4ParticleDefinition* type = 0; 518 G4ParticleDefinition* type = 0;
519 519
520 // Fluorescence data start from element 6 520 // Fluorescence data start from element 6
521 521
522 if (Fluorescence() && Z > 5 && (bindingEnerg 522 if (Fluorescence() && Z > 5 && (bindingEnergy >= cutForPhotons
523 || bindingEnergy >= cutForElectro 523 || bindingEnergy >= cutForElectrons)) {
524 524
525 secondaryVector = deexcitationManager.Gene 525 secondaryVector = deexcitationManager.GenerateParticles(Z, shellId);
526 526
527 if (secondaryVector != 0) { 527 if (secondaryVector != 0) {
528 528
529 nSecondaries = secondaryVector->size(); 529 nSecondaries = secondaryVector->size();
530 for (size_t i = 0; i<nSecondaries; i++) 530 for (size_t i = 0; i<nSecondaries; i++) {
531 531
532 aSecondary = (*secondaryVector)[i]; 532 aSecondary = (*secondaryVector)[i];
533 if (aSecondary) { 533 if (aSecondary) {
534 534
535 G4double e = aSecondary->GetKineticE 535 G4double e = aSecondary->GetKineticEnergy();
536 type = aSecondary->GetDefinition(); 536 type = aSecondary->GetDefinition();
537 if (e < theEnergyDeposit && 537 if (e < theEnergyDeposit &&
538 ((type == G4Gamma::Gamma() && 538 ((type == G4Gamma::Gamma() && e > cutForPhotons ) ||
539 (type == G4Electron::Electron 539 (type == G4Electron::Electron() && e > cutForElectrons ))) {
540 540
541 theEnergyDeposit -= e; 541 theEnergyDeposit -= e;
542 totalNumber++; 542 totalNumber++;
543 543
544 } else { 544 } else {
545 545
546 delete aSecondary; 546 delete aSecondary;
547 (*secondaryVector)[i] = 0; 547 (*secondaryVector)[i] = 0;
548 } 548 }
549 } 549 }
550 } 550 }
551 } 551 }
552 } 552 }
553 553
554 // Save delta-electrons 554 // Save delta-electrons
555 555
556 aParticleChange.SetNumberOfSecondaries(total 556 aParticleChange.SetNumberOfSecondaries(totalNumber);
557 aParticleChange.AddSecondary(theDeltaRay); 557 aParticleChange.AddSecondary(theDeltaRay);
558 558
559 // Save Fluorescence and Auger 559 // Save Fluorescence and Auger
560 560
561 if (secondaryVector) { 561 if (secondaryVector) {
562 562
563 for (size_t l = 0; l < nSecondaries; l++) 563 for (size_t l = 0; l < nSecondaries; l++) {
564 564
565 aSecondary = (*secondaryVector)[l]; 565 aSecondary = (*secondaryVector)[l];
566 566
567 if(aSecondary) { 567 if(aSecondary) {
568 568
569 aParticleChange.AddSecondary(aSecondar 569 aParticleChange.AddSecondary(aSecondary);
570 } 570 }
571 } 571 }
572 delete secondaryVector; 572 delete secondaryVector;
573 } 573 }
574 574
575 if(theEnergyDeposit < 0.) { 575 if(theEnergyDeposit < 0.) {
576 G4cout << "G4LowEnergyIonisation: Negative 576 G4cout << "G4LowEnergyIonisation: Negative energy deposit: "
577 << theEnergyDeposit/eV << " eV" << 577 << theEnergyDeposit/eV << " eV" << G4endl;
578 theEnergyDeposit = 0.0; 578 theEnergyDeposit = 0.0;
579 } 579 }
580 aParticleChange.ProposeLocalEnergyDeposit(th 580 aParticleChange.ProposeLocalEnergyDeposit(theEnergyDeposit);
581 581
582 return G4VContinuousDiscreteProcess::PostSte 582 return G4VContinuousDiscreteProcess::PostStepDoIt(track, step);
583 } 583 }
584 584
585 585
586 void G4LowEnergyIonisation::PrintInfoDefinitio 586 void G4LowEnergyIonisation::PrintInfoDefinition()
587 { 587 {
588 G4String comments = "Total cross sections fr 588 G4String comments = "Total cross sections from EEDL database.";
589 comments += "\n Gamma energy sampled fr 589 comments += "\n Gamma energy sampled from a parametrised formula.";
590 comments += "\n Implementation of the c 590 comments += "\n Implementation of the continuous dE/dx part.";
591 comments += "\n At present it can be us 591 comments += "\n At present it can be used for electrons ";
592 comments += "in the energy range [250eV,100G 592 comments += "in the energy range [250eV,100GeV].";
593 comments += "\n The process must work w 593 comments += "\n The process must work with G4LowEnergyBremsstrahlung.";
594 594
595 G4cout << G4endl << GetProcessName() << ": 595 G4cout << G4endl << GetProcessName() << ": " << comments << G4endl;
596 } 596 }
597 597
598 G4bool G4LowEnergyIonisation::IsApplicable(con 598 G4bool G4LowEnergyIonisation::IsApplicable(const G4ParticleDefinition& particle)
599 { 599 {
600 return ( (&particle == G4Electron::Electron 600 return ( (&particle == G4Electron::Electron()) );
601 } 601 }
602 602
603 std::vector<G4DynamicParticle*>* 603 std::vector<G4DynamicParticle*>*
604 G4LowEnergyIonisation::DeexciteAtom(const G4Ma 604 G4LowEnergyIonisation::DeexciteAtom(const G4MaterialCutsCouple* couple,
605 G4double incidentEnerg 605 G4double incidentEnergy,
606 G4double eLoss) 606 G4double eLoss)
607 { 607 {
608 // create vector of secondary particles 608 // create vector of secondary particles
609 const G4Material* material = couple->GetMate 609 const G4Material* material = couple->GetMaterial();
610 610
611 std::vector<G4DynamicParticle*>* partVector 611 std::vector<G4DynamicParticle*>* partVector =
612 new std::vect 612 new std::vector<G4DynamicParticle*>;
613 613
614 if(eLoss > cutForPhotons && eLoss > cutForEl 614 if(eLoss > cutForPhotons && eLoss > cutForElectrons) {
615 615
616 const G4RDAtomicTransitionManager* transit 616 const G4RDAtomicTransitionManager* transitionManager =
617 G4RDAtomicTrans 617 G4RDAtomicTransitionManager::Instance();
618 618
619 size_t nElements = material->GetNumberOfEl 619 size_t nElements = material->GetNumberOfElements();
620 const G4ElementVector* theElementVector = 620 const G4ElementVector* theElementVector = material->GetElementVector();
621 621
622 std::vector<G4DynamicParticle*>* secVector 622 std::vector<G4DynamicParticle*>* secVector = 0;
623 G4DynamicParticle* aSecondary = 0; 623 G4DynamicParticle* aSecondary = 0;
624 G4ParticleDefinition* type = 0; 624 G4ParticleDefinition* type = 0;
625 G4double e; 625 G4double e;
626 G4ThreeVector position; 626 G4ThreeVector position;
627 G4int shell, shellId; 627 G4int shell, shellId;
628 628
629 // sample secondaries 629 // sample secondaries
630 630
631 G4double eTot = 0.0; 631 G4double eTot = 0.0;
632 std::vector<G4int> n = 632 std::vector<G4int> n =
633 shellVacancy->GenerateNumberOfIonis 633 shellVacancy->GenerateNumberOfIonisations(couple,
634 634 incidentEnergy,eLoss);
635 for (size_t i=0; i<nElements; i++) { 635 for (size_t i=0; i<nElements; i++) {
636 636
637 G4int Z = (G4int)((*theElementVector)[i] 637 G4int Z = (G4int)((*theElementVector)[i]->GetZ());
638 size_t nVacancies = n[i]; 638 size_t nVacancies = n[i];
639 639
640 G4double maxE = transitionManager->Shell 640 G4double maxE = transitionManager->Shell(Z, 0)->BindingEnergy();
641 641
642 if (nVacancies && Z > 5 && (maxE>cutForP 642 if (nVacancies && Z > 5 && (maxE>cutForPhotons || maxE>cutForElectrons)) {
643 643
644 for (size_t j=0; j<nVacancies; j++) { 644 for (size_t j=0; j<nVacancies; j++) {
645 645
646 shell = crossSectionHandler->SelectRandomS 646 shell = crossSectionHandler->SelectRandomShell(Z, incidentEnergy);
647 shellId = transitionManager->Shell(Z 647 shellId = transitionManager->Shell(Z, shell)->ShellId();
648 G4double maxEShell = 648 G4double maxEShell =
649 transitionManager->Shell( 649 transitionManager->Shell(Z, shell)->BindingEnergy();
650 650
651 if (maxEShell>cutForPhotons || maxES 651 if (maxEShell>cutForPhotons || maxEShell>cutForElectrons ) {
652 652
653 secVector = deexcitationManager.Generate 653 secVector = deexcitationManager.GenerateParticles(Z, shellId);
654 654
655 if (secVector != 0) { 655 if (secVector != 0) {
656 656
657 for (size_t l = 0; l<secVector->size() 657 for (size_t l = 0; l<secVector->size(); l++) {
658 658
659 aSecondary = (*secVector)[l]; 659 aSecondary = (*secVector)[l];
660 if (aSecondary != 0) { 660 if (aSecondary != 0) {
661 661
662 e = aSecondary->GetKineticEnergy() 662 e = aSecondary->GetKineticEnergy();
663 type = aSecondary->GetDefinition() 663 type = aSecondary->GetDefinition();
664 if ( eTot + e <= eLoss && 664 if ( eTot + e <= eLoss &&
665 ((type == G4Gamma::Gamma() && e 665 ((type == G4Gamma::Gamma() && e>cutForPhotons ) ||
666 (type == G4Electron::Electron() 666 (type == G4Electron::Electron() && e>cutForElectrons))) {
667 667
668 eTot += e; 668 eTot += e;
669 partVector->push_bac 669 partVector->push_back(aSecondary);
670 670
671 } else { 671 } else {
672 672
673 delete aSecondary; 673 delete aSecondary;
674 674
675 } 675 }
676 } 676 }
677 } 677 }
678 delete secVector; 678 delete secVector;
679 } 679 }
680 } 680 }
681 } 681 }
682 } 682 }
683 } 683 }
684 } 684 }
685 return partVector; 685 return partVector;
686 } 686 }
687 687
688 G4double G4LowEnergyIonisation::GetMeanFreePat 688 G4double G4LowEnergyIonisation::GetMeanFreePath(const G4Track& track,
689 G4double , // previousStepSize 689 G4double , // previousStepSize
690 G4ForceCondition* cond) 690 G4ForceCondition* cond)
691 { 691 {
692 *cond = NotForced; 692 *cond = NotForced;
693 G4int index = (track.GetMaterialCutsCouple( 693 G4int index = (track.GetMaterialCutsCouple())->GetIndex();
694 const G4RDVEMDataSet* data = theMeanFreePat 694 const G4RDVEMDataSet* data = theMeanFreePath->GetComponent(index);
695 G4double meanFreePath = data->FindValue(tra 695 G4double meanFreePath = data->FindValue(track.GetKineticEnergy());
696 return meanFreePath; 696 return meanFreePath;
697 } 697 }
698 698
699 void G4LowEnergyIonisation::SetCutForLowEnSecP 699 void G4LowEnergyIonisation::SetCutForLowEnSecPhotons(G4double cut)
700 { 700 {
701 cutForPhotons = cut; 701 cutForPhotons = cut;
702 deexcitationManager.SetCutForSecondaryPhoton 702 deexcitationManager.SetCutForSecondaryPhotons(cut);
703 } 703 }
704 704
705 void G4LowEnergyIonisation::SetCutForLowEnSecE 705 void G4LowEnergyIonisation::SetCutForLowEnSecElectrons(G4double cut)
706 { 706 {
707 cutForElectrons = cut; 707 cutForElectrons = cut;
708 deexcitationManager.SetCutForAugerElectrons( 708 deexcitationManager.SetCutForAugerElectrons(cut);
709 } 709 }
710 710
711 void G4LowEnergyIonisation::ActivateAuger(G4bo 711 void G4LowEnergyIonisation::ActivateAuger(G4bool val)
712 { 712 {
713 deexcitationManager.ActivateAugerElectronPro 713 deexcitationManager.ActivateAugerElectronProduction(val);
714 } 714 }
715 715
716 716