Context Navigation

G4LowEnergyIonisation.cc @ 1347

Last change on this file since 1347 was 1347, checked in by garnier, 13 years ago
geant4 tag 9.4
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-04-ref-00 $
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
123	G4LowEnergyIonisation::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
146	G4LowEnergyIonisation::~G4LowEnergyIonisation()
147	{
148	delete crossSectionHandler;
149	delete energySpectrum;
150	delete theMeanFreePath;
151	delete shellVacancy;
152	}
153
154
155	void 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
241	void 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
413	G4VParticleChange* 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.0electron_mass_c2));
453	G4double primaryMom = std::sqrt(primaryKinE(primaryKinE + 2.0electron_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.0electron_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 = primaryMomprimaryDirection.x() - deltaMomdirx;
488	G4double finalPy = primaryMomprimaryDirection.y() - deltaMomdiry;
489	G4double finalPz = primaryMomprimaryDirection.z() - deltaMomdirz;
490
491	// create G4DynamicParticle object for delta ray
492	G4DynamicParticle* theDeltaRay = new G4DynamicParticle();
493	theDeltaRay->SetKineticEnergy(tDelta);
494	G4double norm = 1.0/std::sqrt(dirxdirx + dirydiry + 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(finalPxfinalPx+finalPyfinalPy+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
596	void 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
608	G4bool G4LowEnergyIonisation::IsApplicable(const G4ParticleDefinition& particle)
609	{
610	return ( (&particle == G4Electron::Electron()) );
611	}
612
613	std::vector<G4DynamicParticle>
614	G4LowEnergyIonisation::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
698	G4double 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
709	void G4LowEnergyIonisation::SetCutForLowEnSecPhotons(G4double cut)
710	{
711	cutForPhotons = cut;
712	deexcitationManager.SetCutForSecondaryPhotons(cut);
713	}
714
715	void G4LowEnergyIonisation::SetCutForLowEnSecElectrons(G4double cut)
716	{
717	cutForElectrons = cut;
718	deexcitationManager.SetCutForAugerElectrons(cut);
719	}
720
721	void G4LowEnergyIonisation::ActivateAuger(G4bool val)
722	{
723	deexcitationManager.ActivateAugerElectronProduction(val);
724	}
725

Note: See TracBrowser for help on using the repository browser.

Context Navigation

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

Download in other formats: