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

Last change on this file since 1347 was 1347, checked in by garnier, 13 years ago

geant4 tag 9.4
File size: 10.4 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: G4OrlicLiCrossSection.cc,v 1.6 2010/11/22 18:32:00 mantero Exp $
27// GEANT4 tag $Name: geant4-09-04-ref-00 $
28//
29// Author: Haifa Ben Abdelouahed
30//         
31//
32// History:
33// -----------
34//  23 Apr 2008   H. Ben Abdelouahed   1st implementation
35//  28 Apr 2008   MGP        Major revision according to a design iteration
36//  21 Apr 2009   ALF Some correction for compatibility to G4VShellCrossSection
37//                and changed name to G4OrlicLiCrossSection
38//
39// -------------------------------------------------------------------
40// Class description:
41// Low Energy Electromagnetic Physics, Cross section, proton ionisation, L shell
42// Further documentation available from http://www.ge.infn.it/geant4/lowE
43// -------------------------------------------------------------------
44
45
46#include "globals.hh"
47#include "G4OrlicLiCrossSection.hh"
48#include "G4Proton.hh"
49
50
51G4OrlicLiCrossSection::G4OrlicLiCrossSection()
52{ 
53
54  transitionManager =  G4AtomicTransitionManager::Instance();
55
56}
57
58G4OrlicLiCrossSection::~G4OrlicLiCrossSection()
59{ 
60
61}
62
63//this L-CrossSection calculation method is done according to
64//I.ORLIC, C.H.SOW and S.M.TANG,International Journal of PIXE.Vol.4(1997) 217-230       
65 
66
67//*****************************************************************************************************************************************
68
69G4double G4OrlicLiCrossSection::CalculateL1CrossSection(G4int zTarget, G4double energyIncident)
70                               
71{
72
73  if ( (energyIncident < 0.1*MeV) || energyIncident > 10*MeV )
74
75    {return 0;}
76
77
78
79  G4double  massIncident; 
80
81  G4Proton* aProtone = G4Proton::Proton();
82   
83   massIncident = aProtone->GetPDGMass(); 
84 
85  G4double l1BindingEnergy = (transitionManager->Shell(zTarget,1)->BindingEnergy())/keV;
86 
87  G4double lamda =  massIncident/electron_mass_c2;
88
89  G4double normalizedEnergy =  (energyIncident/keV)/(lamda*l1BindingEnergy);
90
91  G4double x = std::log(normalizedEnergy);
92
93  G4double a0 = 0.;
94  G4double a1 = 0.;
95  G4double a2 = 0.;
96  G4double a3 = 0.;
97  G4double a4 = 0.;
98  G4double a5 = 0.;
99  G4double a6 = 0.; 
100  G4double a7 = 0.; 
101  G4double a8 = 0.; 
102  G4double a9 = 0.;
103
104                                       
105   if ( zTarget>=14 && zTarget<=40) 
106    {
107
108      return 0;
109      /*
110      // parameters used for calculating total L cross section
111      a0=12.5081;
112      a1=0.2177;
113      a2=-0.3758;
114      a3=0.0096;
115      a4=0.0073;
116      a5=0.0022;
117      a6=0.;
118      a7=0.;
119      a8=0.;
120      a9=0.; */
121    }
122  else 
123    { 
124       if ( zTarget>=41 &&  zTarget<=50 )
125       {         
126      a0=11.274881;
127      a1=-0.187401;
128      a2=-0.943341;
129      a3=-1.47817;
130      a4=-1.282343;
131      a5=-0.386544; 
132      a6=-0.037932;
133      a7=0.;
134      a8=0.;
135      a9=0.;
136        }
137     
138      else 
139        {
140          if ( zTarget>=51 &&  zTarget<=60 ) 
141            { 
142              a0=11.242637;
143              a1=-0.162515;
144              a2=1.035774;
145              a3=3.970908;
146              a4=3.968233;
147              a5=1.655714;
148              a6=0.058885;
149              a7=-0.155743;
150              a8=-0.042228;
151              a9=-0.003371; 
152            } 
153       
154          else 
155            {
156              if ( zTarget>=61 &&  zTarget<=70 ) 
157                { 
158                  a0=6.476722;
159                  a1=-25.804787;
160                  a2=-54.061629;
161                  a3=-56.684589;
162                  a4=-33.223367;
163                  a5=-11.034979; 
164                  a6=-2.042851;
165                  a7=-0.194075;
166                  a8=-0.007252;
167                  a9=0.;
168                } 
169              else 
170                {
171                  if ( zTarget>=71 &&  zTarget<=80 ) 
172                    { 
173                      a0=12.776794;
174                      a1=6.562907;
175                      a2=10.158703;
176                      a3=7.432592;
177                      a4=2.332036;
178                      a5=0.317946; 
179                      a6=0.014479;
180                      a7=0.;
181                      a8=0.;
182                      a9=0.;
183                    } 
184                  else 
185                    {
186                      if ( zTarget>=81 &&  zTarget<=92 ) 
187                        { 
188                          a0=28.243087;
189                          a1=50.199585;
190                          a2=58.281684;
191                          a3=34.130538;
192                          a4=10.268531;
193                          a5=1.525302; 
194                          a6=0.08835;
195                          a7=0.;
196                          a8=0.;
197                          a9=0.;
198                        }
199                      else
200                        { 
201                          G4cout << "ERROR: L1 Cross-Section exist only for ZTarget between 14 and 92!!! " << G4endl;
202                       
203                        }
204                    }
205                }
206            }
207        }
208      }
209     
210
211G4double analyticalFunction = a0 + (a1*x)+(a2*x*x)+(a3*std::pow(x,3))+(a4*std::pow(x,4))+(a5*std::pow(x,5))+(a6*std::pow(x,6))+
212        (a7*std::pow(x,7))+(a8*std::pow(x,8))+(a9*std::pow(x,9)); 
213
214
215
216  G4double L1crossSection =  std::exp(analyticalFunction)/(l1BindingEnergy*l1BindingEnergy); 
217
218
219  if (L1crossSection >= 0) {
220    return L1crossSection * barn;
221  }
222  else {return 0;}
223
224}
225
226//*****************************************************************************************************************************************
227
228
229G4double G4OrlicLiCrossSection::CalculateL2CrossSection(G4int zTarget, G4double energyIncident)
230                               
231{
232
233
234  if ( (energyIncident < 0.1*MeV) || energyIncident > 10*MeV )
235
236    {return 0;}
237
238
239  G4double  massIncident; 
240
241  G4Proton* aProtone = G4Proton::Proton();
242   
243   massIncident = aProtone->GetPDGMass(); 
244 
245 G4double L2crossSection;
246
247 if (zTarget<41)
248   { 
249     L2crossSection =0.;
250   }
251 else
252   {
253
254 G4double l2BindingEnergy = (transitionManager->Shell(zTarget,2)->BindingEnergy())/keV;
255 
256  G4double lamda =  massIncident/electron_mass_c2;
257
258  G4double normalizedEnergy =  (energyIncident/keV)/(lamda*l2BindingEnergy);
259
260  G4double x = std::log(normalizedEnergy);
261
262  G4double a0 = 0.;
263  G4double a1 = 0.;
264  G4double a2 = 0.;
265  G4double a3 = 0.;
266  G4double a4 = 0.;
267  G4double a5 = 0.;
268 
269      if ( zTarget>=41 &&  zTarget<=50 ) 
270        { 
271      a0=11.194798;
272      a1=0.178807;
273      a2=-0.449865;
274      a3=-0.063528;
275      a4=-0.015364;
276      a5=0.; 
277        }
278     
279      else 
280        {
281          if ( zTarget>=51 &&  zTarget<=60 ) 
282            { 
283              a0=11.241409;
284              a1=0.149635;
285              a2=-0.633269;
286              a3=-0.17834;
287              a4=-0.034743;
288              a5=0.006474;
289
290            } 
291       
292          else 
293            {
294              if ( zTarget>=61 &&  zTarget<=70 ) 
295                { 
296                  a0=11.247424;
297                  a1=0.203051;
298                  a2=-0.219083;
299                  a3=0.164514;
300                  a4=0.058692;
301                  a5=0.007866; 
302                } 
303              else 
304                {
305                  if ( zTarget>=71 &&  zTarget<=80 ) 
306                    { 
307                      a0=11.229924;
308                      a1=-0.087241;
309                      a2=-0.753908;
310                      a3=-0.181546;
311                      a4=-0.030406;
312                      a5=0.; 
313                    } 
314                  else 
315                    {
316                      if ( zTarget>=81 &&  zTarget<=92 ) 
317                        { 
318                          a0=11.586671;
319                          a1=0.730838;
320                          a2=-0.056713;
321                          a3=0.053262;
322                          a4=-0.003672;
323                          a5=0.; 
324                        }
325                      else
326                        { 
327                          G4cout << "ERROR: L2 Cross-Section exist only for ZTarget between 14 and 92!!! " << G4endl;
328                       
329                        }
330                    }
331                }
332            }
333        }
334     
335
336 G4double analyticalFunction = a0 + (a1*x)+(a2*x*x)+(a3*std::pow(x,3))+(a4*std::pow(x,4))+(a5*std::pow(x,5)); 
337
338
339   L2crossSection =  std::exp(analyticalFunction)/(l2BindingEnergy*l2BindingEnergy); 
340                                                                           
341   } 
342
343  if (L2crossSection >= 0) {
344    return L2crossSection * barn;
345  }
346  else {return 0;}
347   
348}
349
350//*****************************************************************************************************************************************
351
352
353G4double G4OrlicLiCrossSection::CalculateL3CrossSection(G4int zTarget, G4double energyIncident)
354                               
355{
356
357  if ( (energyIncident < 0.1*MeV) || energyIncident > 10*MeV )
358
359    {return 0;}
360
361
362
363  G4double  massIncident; 
364
365  G4Proton* aProtone = G4Proton::Proton();
366   
367  massIncident = aProtone->GetPDGMass(); 
368
369 
370 G4double L3crossSection;
371
372 if (zTarget<41)
373   { 
374     L3crossSection =0.;
375   }
376 else
377   {
378
379 G4double l3BindingEnergy = (transitionManager->Shell(zTarget,3)->BindingEnergy())/keV;
380 
381 
382  G4double lamda =  massIncident/electron_mass_c2;
383
384  G4double normalizedEnergy =  (energyIncident/keV)/(lamda*l3BindingEnergy);
385
386  G4double x = std::log(normalizedEnergy);
387
388
389  G4double a0 = 0.;
390  G4double a1 = 0.;
391  G4double a2 = 0.;
392  G4double a3 = 0.;
393  G4double a4 = 0.;
394 
395      if ( zTarget>=41 &&  zTarget<=50 ) 
396        { 
397      a0=11.91837;
398      a1=0.03064;
399      a2=-0.657644;
400      a3=-0.14532;
401      a4=-0.026059;     
402        }
403     
404      else 
405        {
406          if ( zTarget>=51 &&  zTarget<=60 ) 
407            { 
408              a0=11.909485;
409              a1=0.15918;
410              a2=-0.588004;
411              a3=-0.159466;
412              a4=-0.033184;           
413            } 
414       
415          else 
416            {
417              if ( zTarget>=61 &&  zTarget<=70 ) 
418                { 
419                  a0=11.878472;
420                  a1=-0.137007;
421                  a2=-0.959475;
422                  a3=-0.316505;
423                  a4=-0.054154; 
424                } 
425              else 
426                {
427                  if ( zTarget>=71 &&  zTarget<=80 ) 
428                    { 
429                      a0=11.802538;
430                      a1=-0.371796;
431                      a2=-1.052238;
432                      a3=-0.28766;
433                      a4=-0.042608;         
434                    } 
435                  else 
436                    {
437                      if ( zTarget>=81 &&  zTarget<=92 ) 
438                        { 
439                          a0=11.423712;
440                          a1=-1.428823;
441                          a2=-1.946979;
442                          a3=-0.585198;
443                          a4=-0.076467;                 
444                        }
445                      else
446                        { 
447                          G4cout << "ERROR: L3 Cross-Section exist only for ZTarget between 14 and 92!!! " << G4endl;
448                       
449                        }
450                    }
451                }
452            }
453        }
454     
455
456 G4double analyticalFunction = a0 + (a1*x)+(a2*x*x)+(a3*std::pow(x,3))+(a4*std::pow(x,4)); 
457
458
459   L3crossSection =  std::exp(analyticalFunction)/(l3BindingEnergy*l3BindingEnergy); 
460                                                                           
461   } 
462  if (L3crossSection >= 0) {
463    return L3crossSection * barn;
464  }
465  else {return 0;}
466
467
468}
Note: See TracBrowser for help on using the repository browser.