source: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionIonisationRudd.cc@ 830

Last change on this file since 830 was 819, checked in by garnier, 17 years ago

import all except CVS
  • Property svn:executable set to *
File size: 11.4 KB
RevLine 
[819]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//
27// $Id: G4CrossSectionIonisationRudd.cc,v 1.2 2007/11/09 16:20:16 pia Exp $
28// GEANT4 tag $Name: $
29//
30// Contact Author: Maria Grazia Pia (Maria.Grazia.Pia@cern.ch)
31//
32// Reference: TNS Geant4-DNA paper
33// Reference for implementation model: NIM. 155, pp. 145-156, 1978
34
35// History:
36// -----------
37// Date Name Modification
38// 28 Apr 2007 M.G. Pia Created in compliance with design described in TNS paper
39//
40// -------------------------------------------------------------------
41
42// Class description:
43// Geant4-DNA Cross total cross section for electron elastic scattering in water
44// Reference: TNS Geant4-DNA paper
45// S. Chauvie et al., Geant4 physics processes for microdosimetry simulation:
46// design foundation and implementation of the first set of models,
47// IEEE Trans. Nucl. Sci., vol. 54, no. 6, Dec. 2007.
48// Further documentation available from http://www.ge.infn.it/geant4/dna
49
50// -------------------------------------------------------------------
51
52
53#include "G4CrossSectionIonisationRudd.hh"
54#include "G4ParticleDefinition.hh"
55#include "G4Electron.hh"
56#include "G4Proton.hh"
57#include "G4Track.hh"
58#include "G4LogLogInterpolation.hh"
59#include "G4SystemOfUnits.hh"
60#include "G4DNAGenericIonsManager.hh"
61
62G4CrossSectionIonisationRudd::G4CrossSectionIonisationRudd()
63{
64
65 name = "IonisationRudd";
66
67 // Default energy limits (defined for protection against anomalous behaviour only)
68 // ZERO LOW ENERGY LIMIT FOR ALLOWED KILLING
69 lowEnergyLimitDefault = 0 * eV;
70 highEnergyLimitDefault = 100 * MeV;
71
72 G4String fileProton("dna/sigma_ionisation_p_rudd");
73 G4String fileHydrogen("dna/sigma_ionisation_h_rudd");
74 G4String fileAlphaPlusPlus("dna/sigma_ionisation_alphaplusplus_rudd");
75 G4String fileAlphaPlus("dna/sigma_ionisation_alphaplus_rudd");
76 G4String fileHelium("dna/sigma_ionisation_he_rudd");
77
78 G4DNAGenericIonsManager *instance;
79 instance = G4DNAGenericIonsManager::Instance();
80 G4ParticleDefinition* protonDef = G4Proton::ProtonDefinition();
81 G4ParticleDefinition* hydrogenDef = instance->GetIon("hydrogen");
82 G4ParticleDefinition* alphaPlusPlusDef = instance->GetIon("alpha++");
83 G4ParticleDefinition* alphaPlusDef = instance->GetIon("alpha+");
84 G4ParticleDefinition* heliumDef = instance->GetIon("helium");
85
86 G4String proton;
87 G4String hydrogen;
88 G4String alphaPlusPlus;
89 G4String alphaPlus;
90 G4String helium;
91
92 // Factor to scale microscopic/macroscopic cross section data in water
93 // ---- MGP ---- Hardcoded (taken from prototype code); to be replaced with proper calculation
94 G4double scaleFactor = 1 * m*m;
95
96 // Data members for protons
97
98 if (protonDef != 0)
99 {
100 proton = protonDef->GetParticleName();
101 tableFile[proton] = fileProton;
102
103 // Energy limits
104 lowEnergyLimit[proton] = 0. * eV;
105 highEnergyLimit[proton] = 500. * keV;
106
107 // Create data set with proton cross section data and load values stored in file
108 G4DNACrossSectionDataSet* tableProton = new G4DNACrossSectionDataSet(new G4LogLogInterpolation,
109 eV,
110 scaleFactor );
111 tableProton->LoadData(fileProton);
112
113 // Insert key-table pair in map
114 tableData[proton] = tableProton;
115 }
116 else
117 {
118 G4Exception("G4CrossSectionIonisationRudd Constructor: proton is not defined");
119 }
120
121 // Data members for hydrogen
122
123 if (hydrogenDef != 0)
124 {
125 hydrogen = hydrogenDef->GetParticleName();
126 tableFile[hydrogen] = fileHydrogen;
127
128 // Energy limits
129 lowEnergyLimit[hydrogen] = 0. * eV;
130 highEnergyLimit[hydrogen] = 100. * MeV;
131
132 // Create data set with hydrogen cross section data and load values stored in file
133 G4DNACrossSectionDataSet* tableHydrogen = new G4DNACrossSectionDataSet(new G4LogLogInterpolation,
134 eV,
135 scaleFactor );
136 tableHydrogen->LoadData(fileHydrogen);
137
138 // Insert key-table pair in map
139 tableData[hydrogen] = tableHydrogen;
140 }
141 else
142 {
143 G4Exception("G4CrossSectionIonisationRudd Constructor: hydrogen is not defined");
144 }
145
146 // Data members for alphaPlusPlus
147
148 if (alphaPlusPlusDef != 0)
149 {
150 alphaPlusPlus = alphaPlusPlusDef->GetParticleName();
151 tableFile[alphaPlusPlus] = fileAlphaPlusPlus;
152
153 // Energy limits
154 lowEnergyLimit[alphaPlusPlus] = 0. * keV;
155 highEnergyLimit[alphaPlusPlus] = 10. * MeV;
156
157 // Create data set with hydrogen cross section data and load values stored in file
158 G4DNACrossSectionDataSet* tableAlphaPlusPlus = new G4DNACrossSectionDataSet(new G4LogLogInterpolation,
159 eV,
160 scaleFactor );
161 tableAlphaPlusPlus->LoadData(fileAlphaPlusPlus);
162
163 // Insert key-table pair in map
164 tableData[alphaPlusPlus] = tableAlphaPlusPlus;
165 }
166 else
167 {
168 G4Exception("G4CrossSectionIonisationRudd Constructor: alphaPlusPlus is not defined");
169 }
170
171 // Data members for alphaPlus
172
173 if (alphaPlusDef != 0)
174 {
175 alphaPlus = alphaPlusDef->GetParticleName();
176 tableFile[alphaPlus] = fileAlphaPlus;
177
178 // Energy limits
179 lowEnergyLimit[alphaPlus] = 0. * keV;
180 highEnergyLimit[alphaPlus] = 10. * MeV;
181
182 // Create data set with hydrogen cross section data and load values stored in file
183 G4DNACrossSectionDataSet* tableAlphaPlus = new G4DNACrossSectionDataSet(new G4LogLogInterpolation,
184 eV,
185 scaleFactor );
186 tableAlphaPlus->LoadData(fileAlphaPlus);
187
188 // Insert key-table pair in map
189 tableData[alphaPlus] = tableAlphaPlus;
190 }
191 else
192 {
193 G4Exception("G4CrossSectionIonisationRudd Constructor: alphaPlus is not defined");
194 }
195
196 // Data members for helium
197
198 if (heliumDef != 0)
199 {
200 helium = heliumDef->GetParticleName();
201 tableFile[helium] = fileHelium;
202
203 // Energy limits
204 lowEnergyLimit[helium] = 0. * keV;
205 highEnergyLimit[helium] = 10. * MeV;
206
207 // Create data set with hydrogen cross section data and load values stored in file
208 G4DNACrossSectionDataSet* tableHelium = new G4DNACrossSectionDataSet(new G4LogLogInterpolation,
209 eV,
210 scaleFactor );
211 tableHelium->LoadData(fileHelium);
212
213 // Insert key-table pair in map
214 tableData[helium] = tableHelium;
215 }
216 else
217 {
218 G4Exception("G4CrossSectionIonisationRudd Constructor: helium is not defined");
219 }
220}
221
222
223G4CrossSectionIonisationRudd::~G4CrossSectionIonisationRudd()
224{
225 // Destroy the content of the map
226 std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos;
227 for (pos = tableData.begin(); pos != tableData.end(); ++pos)
228 {
229 G4DNACrossSectionDataSet* table = pos->second;
230 delete table;
231 }
232}
233
234
235
236G4double G4CrossSectionIonisationRudd::CrossSection(const G4Track& track )
237{
238 const G4DynamicParticle* particle = track.GetDynamicParticle();
239 G4double k = particle->GetKineticEnergy();
240
241 G4DNAGenericIonsManager *instance;
242 instance = G4DNAGenericIonsManager::Instance();
243
244 if (
245 track.GetDefinition() != G4Proton::ProtonDefinition()
246 &&
247 track.GetDefinition() != instance->GetIon("hydrogen")
248 &&
249 track.GetDefinition() != instance->GetIon("alpha++")
250 &&
251 track.GetDefinition() != instance->GetIon("alpha+")
252 &&
253 track.GetDefinition() != instance->GetIon("helium")
254 )
255
256 G4Exception("G4CrossSectionIonisationRudd: attempting to calculate cross section for wrong particle");
257
258 // Cross section = 0 outside the energy validity limits set in the constructor
259 G4double sigma = 0.;
260
261 // ---- MGP ---- Better handling of these limits to be set in a following design iteration
262 G4double lowLim = lowEnergyLimitDefault;
263 G4double highLim = highEnergyLimitDefault;
264
265 const G4String& particleName = particle->GetDefinition()->GetParticleName();
266
267 // Retrieve energy limits for the current particle type
268
269 std::map< G4String,G4double,std::less<G4String> >::iterator pos1;
270 pos1 = lowEnergyLimit.find(particleName);
271
272 // Lower limit
273 if (pos1 != lowEnergyLimit.end())
274 {
275 lowLim = pos1->second;
276 }
277
278 // Upper limit
279 std::map< G4String,G4double,std::less<G4String> >::iterator pos2;
280 pos2 = highEnergyLimit.find(particleName);
281
282 if (pos2 != highEnergyLimit.end())
283 {
284 highLim = pos2->second;
285 }
286
287 // Verify that the current track is within the energy limits of validity of the cross section model
288 if (k >= lowLim && k <= highLim)
289 {
290 std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos;
291 pos = tableData.find(particleName);
292
293 if (pos != tableData.end())
294 {
295 G4DNACrossSectionDataSet* table = pos->second;
296 if (table != 0)
297 {
298 // Cross section
299 sigma = table->FindValue(k);
300
301
302 // BEGIN ELECTRON CORRECTION
303 // add ONE or TWO electron-water excitation for alpha+ and helium
304
305 if ( particle->GetDefinition() == instance->GetIon("alpha+")
306 ||
307 particle->GetDefinition() == instance->GetIon("helium")
308 )
309 {
310
311 G4DNACrossSectionDataSet* electronDataset = new G4DNACrossSectionDataSet
312 (new G4LogLogInterpolation, eV, (1./3.343e22)*m*m);
313
314 electronDataset->LoadData("dna/sigma_ionisation_e_born");
315
316 G4double kElectron = k * 0.511/3728;
317
318 if ( particle->GetDefinition() == instance->GetIon("alpha+") )
319 {
320 G4double tmp1 = table->FindValue(k) + electronDataset->FindValue(kElectron);
321 delete electronDataset;
322 return tmp1;
323 }
324
325 if ( particle->GetDefinition() == instance->GetIon("helium") )
326 {
327 G4double tmp2 = table->FindValue(k) + 2. * electronDataset->FindValue(kElectron);
328 delete electronDataset;
329 return tmp2;
330 }
331 }
332
333 // END ELECTRON CORRECTION
334
335 }
336 }
337 else
338 {
339 // The track does not corresponds to a particle pertinent to this model
340 G4Exception("G4CrossSectionIonisationRudd: attempting to calculate cross section for wrong particle");
341 }
342 }
343
344 return sigma;
345}
346
Note: See TracBrowser for help on using the repository browser.