Bugfix in AliPoints2Memory
[u/mrichter/AliRoot.git] / GEANT321 / fluka / rbkekv.F
1 *
2 * $Id$
3 *
4 * $Log$
5 * Revision 1.1.1.1  1995/10/24 10:19:58  cernlib
6 * Geant
7 *
8 *
9 #include "geant321/pilot.h"
10 *CMZ :  3.21/02 29/03/94  15.41.44  by  S.Giani
11 *-- Author :
12 *$ CREATE RBKEKV.FOR
13 *COPY RBKEKV
14 *                                                                      *
15 *=== rbkekv ===========================================================*
16 *                                                                      *
17       SUBROUTINE RBKEKV ( IT, EXSOP, TO, AMSS, TKIN, TSTRCK,
18      &                    PSTRCK, ARECL, TKRECL, COD, SID )
19  
20 #include "geant321/dblprc.inc"
21 #include "geant321/dimpar.inc"
22 #include "geant321/iounit.inc"
23 *
24 *----------------------------------------------------------------------*
25 *     version by                     Alfredo Ferrari
26 *                                    INFN - Milan
27 *     last change 19 april 93 by     Alfredo Ferrari
28 *                                    INFN - Milan
29 *
30 *     To be called from the high energy production
31 *
32 *     this is a subroutine of fluka to sample intranuclear cascade
33 *     particles: it is based on the old Rakeka from J. Ranft
34 *
35 *     input variables:
36 *        it    = type of the secondary requested; 1=proton, 2=neutron
37 *        bbtar = atomic weight of the medium
38 *
39 *     output variables:
40 *        tkin  = kinetic energy of the secondary in GeV before applying
41 *                the nuclear well (and eventually the Coulomb barreer)
42 *        tstrck= kinetic energy of the secondary in GeV
43 *        pstrck= momentum of the secondary in GeV/c
44 *        sid,cod = sine and cosine of the angle between
45 *                  the directions of the primary
46 *                  and secondary particles
47 *
48 *********************************************************************
49 *
50       PARAMETER ( PI   = PIPIPI )
51       PARAMETER ( PIO2 = PI / 2.D+00 )
52 #include "geant321/nucdat.inc"
53       LOGICAL LINIT, LZEROI
54       DIMENSION EXSOP (2), AMMOLD (2)
55       DIMENSION PFINIT (50,2), TKINIT (50,2), TSINIT (50,2),
56      &          TRINIT (50,2), NINI (2)
57       REAL RNDM(4)
58       SAVE AMMOLD, PFINIT, TKINIT, TSINIT, TRINIT, NINI
59       DATA AMMOLD / 0.9382796D+00, 0.9395731D+00 /
60 *
61       LINIT = .FALSE.
62       IF ( NINI (IT) .GT. 0 ) THEN
63          PKFRMI = PFINIT (NINI(IT),IT)
64          TKIN   = TKINIT (NINI(IT),IT)
65          TSTRCK = TSINIT (NINI(IT),IT)
66          TKRECL = TRINIT (NINI(IT),IT)
67          NINI  (IT) = NINI (IT) - 1
68          EXSOP (IT) = 0.D+00
69          GO TO 3500
70       ELSE
71          GO TO 200
72       END IF
73       ENTRY RBKINI ( IT, LZEROI, EXSOP, TKIN, TSTRCK,
74      &               PSTRCK, ARECL, TKRECL )
75       LINIT = .TRUE.
76       IF ( LZEROI ) THEN
77          NINI (1) = 0
78          NINI (2) = 0
79          RETURN
80       END IF
81   200 CONTINUE
82 *  In this version the low energy component has been suppressed, since
83 *  they are now produced by the evaporation model, A. Ferrari.
84 *  The parameters needed for the sampling have been already set in
85 *  Corrin:
86       EXSOP (IT) = 0.D+00
87 *  Sample the Fermi momentum
88       ESLPFF = ESLOPE (IT)
89       EXUPFF = EXUPNU (IT)
90       EKUPFF = EKUPNU (IT)
91       EXDWFF = EXMNNU (IT)
92       EKDWFF = EKMNNU (IT)
93       ERCLFF = ERCLAV (IT)
94       IRECNT = 0
95 *  +-------------------------------------------------------------------*
96 *  |
97  100  CONTINUE
98 *  |  Sample the Fermi momentum
99          CALL GRNDM(RNDM,4)
100          PKFRMI = PFRMMX (IT) * MAX ( RNDM (1), RNDM (2),
101      &            RNDM (3) )
102          PKFRSQ = PKFRMI * PKFRMI
103          TKIN   = - ESLPFF * LOG ( EXDWFF - RNDM (4) * (
104      &              EXDWFF - EXUPFF ) )
105          TKRECL = 0.5D+00 * PKFRSQ / ( AMUC12 * ARECL ) * ( 1.D+00 -
106      &            0.25D+00 * PKFRSQ / ( ARECL**2 * AMUCSQ ) )
107          TSTRCK = SQRT ( AMNUSQ (IT) + PKFRSQ ) + TKIN - AMNUCL (IT)
108      &          - V0WELL (IT) - TKRECL - EBNDNG (IT)
109 *  |  +----------------------------------------------------------------*
110 *  |  |
111          IF ( TSTRCK .LE. VEFFNU (IT) - V0WELL (IT) ) THEN
112             EXSOP (IT) = EXSOP (IT) + TKIN
113             IRECNT = IRECNT + 1
114 *  |  |  +-------------------------------------------------------------*
115 *  |  |  |
116             IF ( IRECNT .GT. 10 ) THEN
117 *  |  |  |  +----------------------------------------------------------*
118 *  |  |  |  |
119                IF ( TKRECL - ERCLFF .GT. EKUPFF - EKDWFF .AND.
120      &              IRECNT .LT. 20 ) THEN
121                   ERCLFF = ERCLFF + TKRECL - ERCLFF
122                   EKUPFF = EKUPFF + TKRECL - ERCLFF
123                   EKDWFF = EKDWFF + TKRECL - ERCLFF
124                   AHELP  = EXP ( - ( TKRECL - ERCLFF ) / ESLPFF )
125                   EXUPFF = EXUPFF * AHELP
126                   EXDWFF = EXDWFF * AHELP
127 *  |  |  |  |
128 *  |  |  |  +----------------------------------------------------------*
129 *  |  |  |  |
130                ELSE IF ( IRECNT .GT. 15 ) THEN
131                   TKRECL = MAX ( TKRECL, ERCLFF )
132                   ERCLFF = ERCLFF + TKRECL
133                   EKUPFF = EKUPFF + TKRECL - ERCLFF
134                   EKDWFF = EKDWFF + TKRECL - ERCLFF
135                   AHELP  = EXP ( - TKRECL / ESLPFF )
136                   EXUPFF = EXUPFF * AHELP
137                   EXDWFF = EXDWFF * AHELP
138                END IF
139 *  |  |  |  |
140 *  |  |  |  +----------------------------------------------------------*
141             END IF
142 *  |  |  |
143 *  |  |  +-------------------------------------------------------------*
144             GO TO 100
145 *  +-<|--<--<--<--<--< go to resampling
146          END IF
147 *  |  |
148 *  |  +----------------------------------------------------------------*
149 *  +-------------------------------------------------------------------*
150 *  |
151       IF ( LINIT ) THEN
152          NINI (IT) = NINI (IT) + 1
153          PFINIT (NINI(IT),IT) = PKFRMI
154          TKINIT (NINI(IT),IT) = TKIN
155          TSINIT (NINI(IT),IT) = TSTRCK
156          TRINIT (NINI(IT),IT) = TKRECL
157          RETURN
158       END IF
159 *  |
160 *  +-------------------------------------------------------------------*
161  3500 CONTINUE
162 *  | Masses have been updated
163       PSTRCK = SQRT ( TSTRCK * ( TSTRCK + 2.D+00 * AMNUCL (IT) ) )
164 *
165 ********************* Sample the angle ********************************
166 * Polar angle selection
167       ADE=0.090D0*(1.D0+0.081D0*ATO1O3)/TKIN
168       DEX=EXP(- PIO2 * PIO2 / ADE)
169       AN1=(1.D0-DEX)*ADE/2.D0
170       AN2=DEX*PIO2
171       AN=AN1+AN2
172       AN1=AN1/AN
173       CALL GRNDM(RNDM,1)
174       IF(RNDM(1).GT.AN1)  GO TO 3
175     2 CONTINUE
176       CALL GRNDM(RNDM,1)
177       DE=SQRT(-ADE*LOG(1.D0-RNDM(1)*(1.D0-DEX)))
178       IF(DE.GT.PIO2) GO TO 2
179 C     WRITE(LUNOUT,*)IT,TO,AMSS,SQAMSS,T,P,DE
180       COD = COS (DE)
181       SID = SIN (DE)
182       RETURN
183     3 CONTINUE
184       CALL GRNDM(RNDM,1)
185       COD = - RNDM(1)
186       SID = SQRT ( (1.D0 + COD ) * ( 1.D0 - COD ) )
187 C     WRITE(LUNOUT,*)IT,TO,AMSS,SQAMSS,T,P,DE
188       RETURN
189       ENTRY RBKMIN (IT)
190       NINI(IT) = NINI(IT)-1
191       RETURN
192       END