| Geant4 Cross Reference |
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 // Author: Haifa Ben Abdelouahed
27 //
28 //
29 // History:
30 // -----------
31 // 23 Apr 2008 H. Ben Abdelouahed 1st implementation
32 // 28 Apr 2008 MGP Major revision according to a design iteration
33 // 21 Apr 2009 ALF Some correction for compatibility to G4VShellCrossSection
34 // and changed name to G4OrlicLiCrossSection
35 // 21 Mar 2011 ALF some bug fixing (Z checks, )
36 // 29 Oct 2011 ALF Changed name to G4OrlicLiXsModel
37 //
38 // -------------------------------------------------------------------
39 // Class description:
40 // Low Energy Electromagnetic Physics, Cross section, proton ionisation, L shell
41 // Further documentation available from http://www.ge.infn.it/geant4/lowE
42 // -------------------------------------------------------------------
43
44 #include "G4OrlicLiXsModel.hh"
45
46 #include "globals.hh"
47 #include "G4PhysicalConstants.hh"
48 #include "G4SystemOfUnits.hh"
49 #include "G4Proton.hh"
50 #include "G4Exp.hh"
51
52 //********************************************************************************
53
54 G4OrlicLiXsModel::G4OrlicLiXsModel()
55 {
56 transitionManager = G4AtomicTransitionManager::Instance();
57 }
58
59 //********************************************************************************
60
61 G4OrlicLiXsModel::~G4OrlicLiXsModel()
62 {;}
63
64 //this L-CrossSection calculation method is done according to
65 //I.ORLIC, C.H.SOW and S.M.TANG,International Journal of PIXE.Vol.4(1994) 217-230
66
67 //********************************************************************************
68 G4double G4OrlicLiXsModel::CalculateL1CrossSection(G4int zTarget, G4double energyIncident)
69
70 {
71
72 if ( zTarget < 41 )//fixed: no control on z!
73 {
74 return 0;
75 }
76
77 G4double l1BindingEnergy = transitionManager->Shell(zTarget,1)->BindingEnergy()/keV;
78 G4double lamda = 1836.109; //massIncident/electron_mass_c2;
79
80 G4double normalizedEnergy = (energyIncident/keV)/(lamda*l1BindingEnergy);
81
82 G4double x = std::log(normalizedEnergy);
83
84 G4double a0 = 0.;
85 G4double a1 = 0.;
86 G4double a2 = 0.;
87 G4double a3 = 0.;
88 G4double a4 = 0.;
89 G4double a5 = 0.;
90 G4double a6 = 0.;
91 G4double a7 = 0.;
92 G4double a8 = 0.;
93 G4double a9 = 0.;
94
95 if ( (zTarget>=41 && zTarget<=50) && (normalizedEnergy>=0.013 && normalizedEnergy<=1) )
96 {
97 a0=11.274881;
98 a1=-0.187401;
99 a2=-0.943341;
100 a3=-1.47817;
101 a4=-1.282343;
102 a5=-0.386544;
103 a6=-0.037932;
104 a7=0.;
105 a8=0.;
106 a9=0.;
107 }
108 else if ( (zTarget>=51 && zTarget<=60) && (normalizedEnergy>=0.012 && normalizedEnergy<=0.95))
109 {
110 a0=11.242637;
111 a1=-0.162515;
112 a2=1.035774;
113 a3=3.970908;
114 a4=3.968233;
115 a5=1.655714;
116 a6=0.058885;
117 a7=-0.155743;
118 a8=-0.042228;
119 a9=-0.003371;
120 }
121 else if ( (zTarget>=61 && zTarget<=70) && (normalizedEnergy>=0.01 && normalizedEnergy<=0.6) )
122 {
123 a0=6.476722;
124 a1=-25.804787;
125 a2=-54.061629;
126 a3=-56.684589;
127 a4=-33.223367;
128 a5=-11.034979;
129 a6=-2.042851;
130 a7=-0.194075;
131 a8=-0.007252;
132 a9=0.;
133 }
134 else if ( (zTarget>=71 && zTarget<=80) && (normalizedEnergy>=0.01 && normalizedEnergy<=0.45) )
135 {
136 a0=12.776794;
137 a1=6.562907;
138 a2=10.158703;
139 a3=7.432592;
140 a4=2.332036;
141 a5=0.317946;
142 a6=0.014479;
143 a7=0.;
144 a8=0.;
145 a9=0.;
146 }
147 else if ( (zTarget>=81 && zTarget<=92) && (normalizedEnergy>=0.008 && normalizedEnergy<=0.3) )
148 {
149 a0=28.243087;
150 a1=50.199585;
151 a2=58.281684;
152 a3=34.130538;
153 a4=10.268531;
154 a5=1.525302;
155 a6=0.08835;
156 a7=0.;
157 a8=0.;
158 a9=0.;
159 }
160 else {return 0;}
161
162
163 G4double analyticalFunction = a0 + (a1*x)+(a2*x*x)+(a3*std::pow(x,3))+
164 (a4*std::pow(x,4))+(a5*std::pow(x,5))+(a6*std::pow(x,6))+
165 (a7*std::pow(x,7))+(a8*std::pow(x,8))+(a9*std::pow(x,9));
166
167 G4double L1crossSection = G4Exp(analyticalFunction)/(l1BindingEnergy*l1BindingEnergy);
168
169 if (L1crossSection >= 0) {
170 return L1crossSection * barn;
171 }
172 else {return 0;}
173 }
174
175 //************************************************************************************************
176
177 G4double G4OrlicLiXsModel::CalculateL2CrossSection(G4int zTarget, G4double energyIncident)
178
179 {
180 if ( zTarget < 41) //fixed: no control on z!)
181 {
182 return 0;
183 }
184
185 G4double massIncident;
186 G4Proton* aProtone = G4Proton::Proton();
187
188 massIncident = aProtone->GetPDGMass();
189
190 G4double L2crossSection;
191
192 G4double l2BindingEnergy = (transitionManager->Shell(zTarget,2)->BindingEnergy())/keV;
193 G4double lamda = massIncident/electron_mass_c2;
194 G4double normalizedEnergy = (energyIncident/keV)/(lamda*l2BindingEnergy);
195 G4double x = std::log(normalizedEnergy);
196
197 G4double a0 = 0.;
198 G4double a1 = 0.;
199 G4double a2 = 0.;
200 G4double a3 = 0.;
201 G4double a4 = 0.;
202 G4double a5 = 0.;
203
204 if ( (zTarget>=41 && zTarget<=50) && (normalizedEnergy>=0.015 && normalizedEnergy<=1.5))
205 {
206 a0=11.194798;
207 a1=0.178807;
208 a2=-0.449865;
209 a3=-0.063528;
210 a4=-0.015364;
211 a5=0.;
212 }
213 else if ( (zTarget>=51 && zTarget<=60) && (normalizedEnergy>=0.012 && normalizedEnergy<=1.0))
214 {
215 a0=11.241409;
216 a1=0.149635;
217 a2=-0.633269;
218 a3=-0.17834;
219 a4=-0.034743;
220 a5=0.006474; // a little bit better if this is zero (respect to ecpssr)
221 }
222 else if ( (zTarget>=61 && zTarget<=70) && (normalizedEnergy>=0.01 && normalizedEnergy<=0.65))
223 {
224 a0=11.247424;
225 a1=0.203051;
226 a2=-0.219083;
227 a3=0.164514;
228 a4=0.058692;
229 a5=0.007866;
230 }
231 else if ( (zTarget>=71 && zTarget<=80) && (normalizedEnergy>=0.01 && normalizedEnergy<=0.47))
232 {
233 a0=11.229924;
234 a1=-0.087241;
235 a2=-0.753908;
236 a3=-0.181546;
237 a4=-0.030406;
238 a5=0.;
239 }
240 else if ( (zTarget>=81 && zTarget<=92) && (normalizedEnergy>=0.01 && normalizedEnergy<=0.35))
241 {
242 a0=11.586671;
243 a1=0.730838;
244 a2=-0.056713;
245 a3=0.053262;
246 a4=-0.003672;
247 a5=0.;
248 }
249 else {return 0;}
250
251 G4double analyticalFunction = a0 + (a1*x)+(a2*x*x)+(a3*std::pow(x,3))+
252 (a4*std::pow(x,4))+(a5*std::pow(x,5));
253
254 L2crossSection = G4Exp(analyticalFunction)/(l2BindingEnergy*l2BindingEnergy);
255
256 if (L2crossSection >= 0) {
257 return L2crossSection * barn;
258 }
259 else {return 0;}
260 }
261
262 //**************************************************************************************************
263
264 G4double G4OrlicLiXsModel::CalculateL3CrossSection(G4int zTarget, G4double energyIncident)
265
266 {
267
268 if ( zTarget < 41) //fixed: no control on z!
269 {
270 return 0;
271 }
272
273 G4double massIncident;
274 G4Proton* aProtone = G4Proton::Proton();
275 massIncident = aProtone->GetPDGMass();
276
277 G4double L3crossSection;
278 G4double l3BindingEnergy = (transitionManager->Shell(zTarget,3)->BindingEnergy())/keV;
279
280 G4double lamda = massIncident/electron_mass_c2;
281 G4double normalizedEnergy = (energyIncident/keV)/(lamda*l3BindingEnergy);
282 G4double x = std::log(normalizedEnergy);
283
284 G4double a0 = 0.;
285 G4double a1 = 0.;
286 G4double a2 = 0.;
287 G4double a3 = 0.;
288 G4double a4 = 0.;
289 G4double a5 = 0.;
290
291 if ( (zTarget>=41 && zTarget<=50 ) && (normalizedEnergy>=0.015 && normalizedEnergy<=1.5))
292 {
293 a0=11.91837;
294 a1=0.03064;
295 a2=-0.657644;
296 a3=-0.14532;
297 a4=-0.026059;
298 //a5=-0.044735; Correction to Orlic model as explained in
299 //Abdelhwahed H Incerti S and Mantero A 2009 Nucl. Instrum. Meth.B 267 37
300 }
301 else if ( (zTarget>=51 && zTarget<=60 ) && (normalizedEnergy>=0.013 && normalizedEnergy<=1.1))
302 {
303 a0=11.909485;
304 a1=0.15918;
305 a2=-0.588004;
306 a3=-0.159466;
307 a4=-0.033184;
308 }
309 else if ( (zTarget>=61 && zTarget<=70 ) && (normalizedEnergy>=0.01 && normalizedEnergy<=0.67))
310 {
311 a0=11.878472;
312 a1=-0.137007;
313 a2=-0.959475;
314 a3=-0.316505;
315 a4=-0.054154;
316 }
317 else if ( (zTarget>=71 && zTarget<=80 ) && (normalizedEnergy>=0.013 && normalizedEnergy<=0.5))
318 {
319 a0=11.802538;
320 a1=-0.371796;
321 a2=-1.052238;
322 a3=-0.28766;
323 a4=-0.042608;
324 }
325 else if ( (zTarget>=81 && zTarget<=92 ) && (normalizedEnergy>=0.01 && normalizedEnergy<=0.35))
326 {
327 a0=11.423712;
328 a1=-1.428823;
329 a2=-1.946979;
330 a3=-0.585198;
331 a4=-0.076467;
332 }
333 else {return 0;}
334
335 G4double analyticalFunction = a0 + (a1*x)+(a2*x*x)+
336 (a3*std::pow(x,3))+(a4*std::pow(x,4))+(a5*std::pow(x,5));
337
338 L3crossSection = G4Exp(analyticalFunction)/(l3BindingEnergy*l3BindingEnergy);
339
340 if (L3crossSection >= 0) {
341 return L3crossSection * barn;
342 }
343 else {return 0;}
344 }
345