source: trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyIonisation.cc@ 1201

Last change on this file since 1201 was 1196, checked in by garnier, 16 years ago

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