Bugfix in AliPoints2Memory
[u/mrichter/AliRoot.git] / GEANT321 / gstrag / gkokri.F
CommitLineData
fe4da5cc 1*
2* $Id$
3*
4* $Log$
5* Revision 1.1.1.1 1995/10/24 10:21:37 cernlib
6* Geant
7*
8*
9#include "geant321/pilot.h"
10*CMZ : 3.21/02 29/03/94 15.41.25 by S.Giani
11*-- Author :
12 FUNCTION GKOKRI(E,EMINEV,EMAXEV)
13C.
14C. ******************************************************************
15C. * *
16C. * Input energy in eV, Sandia tables in keV *
17C. * *
18C. * *
19C. ******************************************************************
20C.
21#include "geant321/gcbank.inc"
22#include "geant321/gcjloc.inc"
23#include "geant321/gconsp.inc"
24#include "geant321/gcmate.inc"
25#include "geant321/gc10ev.inc"
26#if !defined(CERNLIB_SINGLE)
27 DOUBLE PRECISION EINV,ECUR,ONE,RES,ZERO,EBEG,EEND
28 DOUBLE PRECISION ALPH,BETA,GAMM,E1,E2,GPSCIN,POLE
29 DOUBLE PRECISION EMAX,EMIN,REST
30#endif
31 PARAMETER (ONE=1,ZERO=0)
32C.
33C. ------------------------------------------------------------------
34C.
35 REST = ZERO
36 EMAX = EMAXEV*1E-3
37 EMIN = EMINEV*1E-3
38 ECUR = E*1E-3
39 EINV = ONE/ECUR
40C Use Sandia data
41 JPHXS = LQ(JPHOT-1)
42 NZ = Q(JPHXS+1)
43 JWEIGH = JPHXS+1+2*NZ
44 DO 30 JZ=1,NZ
45 RES = ZERO
46 EBEG = EMIN
47 WEIGHT = Q(JWEIGH+JZ)
48 JPHFN = LQ(JPHXS-JZ)
49 IPOINT = JPHFN+1
50 IMAX = Q(IPOINT)
51 IPOINT = IPOINT-4
52 DO 10 I = 1,IMAX
53 IPOINT = IPOINT+5
54 EEND = Q(IPOINT)
55 IF(EEND.GT.EMIN) THEN
56 E1 = MAX(EBEG,EMIN)
57 E2 = MIN(EEND,EMAX)
58 J = IPOINT+1
59 IF(ECUR.LE.E2.AND.ECUR.GE.E1) THEN
60*
61* *** The pole of the integration is in this interval
62 EPS1 = (ECUR-EBEG)/ECUR
63 EPS2 = (EEND-ECUR)/EEND
64 IF(EPS1.LT.EPS2) THEN
65*
66* *** First the pole and then a simple integration
67 ALPH = ONE-EPS1
68 E1 = ECUR/ALPH
69 E2 = EEND
70 ELSE
71*
72* *** First a simple integration and then the pole
73 ALPH = ONE-EPS2
74 E1 = EBEG
75 E2 = ECUR*ALPH
76 ENDIF
77 BETA = -LOG(ALPH)
78 GAMM = ONE/ALPH
79 POLE = EINV*(Q(J)*BETA+EINV*( Q(J+1)*(GAMM-ALPH)+
80 + EINV*( Q(J+2)*(0.5*GAMM**2+BETA-0.5*ALPH**2)+EINV*
81 + Q(J+3)*(GAMM**3/3+GAMM-ALPH-ALPH**3/3))))
82 RES = RES - POLE
83 ENDIF
84*
85* *** This is a normal integration
86 RES = RES + GPSCIN(E1,E2,ECUR,Q(J))
87 ENDIF
88 EBEG = EEND
89 IF(EBEG.GE.EMAX) GO TO 20
90 10 CONTINUE
91 20 REST = REST+WEIGHT*RES
92 30 CONTINUE
93C RES value is in cm**2/(g keV)
94 GKOKRI = REST*1E-6*A/AVO
95C Now in Megabarns/eV
96C
97 END