This commit was generated by cvs2svn to compensate for changes in r2,
[u/mrichter/AliRoot.git] / GEANT321 / fluka / fkimpu.F
1 *
2 * $Id$
3 *
4 * $Log$
5 * Revision 1.1.1.1  1995/10/24 10:20:06  cernlib
6 * Geant
7 *
8 *
9 #include "geant321/pilot.h"
10 *CMZ :  3.21/02 29/03/94  15.41.45  by  S.Giani
11 *-- Author :
12 *=== impuls ===========================================================*
13 *
14       SUBROUTINE FKIMPU(HE,HMA,HPS,HPX,HPY,HPZ,LT,LL,B3)
15  
16 #include "geant321/dblprc.inc"
17 #include "geant321/dimpar.inc"
18 #include "geant321/iounit.inc"
19       REAL RNDM(2)
20 *
21 *----------------------------------------------------------------------*
22 *  Impuls89: slight revision by A. Ferrari                             *
23 *----------------------------------------------------------------------*
24 *
25 C*****CHOICE OF THE TRANSVERSAL MOMENTUM AND CALCULATION OF THE
26 C*****LONGITUDINAL MOMENTUM
27 C*****HPS,HPZ ARE THE TRANSVERSAL AND LONGITUDINAL MOMENTUM
28  
29       ESMAX = HE - HMA
30     1 CONTINUE
31 *
32 *  He is the total energy, hma the mass one (input) hpx, hpy, hpz
33 *  the momentum components (output values), hps the transversal
34 *  momentum (output)
35 *
36 *  +-------------------------------------------------------------------*
37 *  |                            cut off for the longitudinal momentum
38       IF (ESMAX .LE. 0.05D0) THEN
39 *  |
40 *  | only transversal momentum!
41 *  |
42          HPS = SQRT(HE**2 - HMA**2)
43          HPZ = 0.D0
44       ELSE
45 *  |
46 *  | hps is the transversal momentum
47 *  |
48          CALL GRNDM(RNDM,2)
49          X = RNDM(1)
50          Y = RNDM(2)
51 *  |
52 *  |  Es is the transverse kinetic energy
53 *  |
54          ES = -2.D0/(B3**2)*LOG(X*Y)
55          IF (ES .GT. ESMAX) GO TO 1
56          HPS = ES**2 + 2.D0*ES*HMA
57          IF (LL .EQ. 0) Y1 = 1.D0
58          IF (LL .EQ. 1) Y1 =-1.D0
59          HPZ = SQRT(HE**2 - HMA**2 - HPS) * Y1
60          HPS = SQRT (HPS)
61       END IF
62 *  |
63 *  +-------------------------------------------------------------------*
64       CALL SFECFE(SFE,CFE)
65       SIP = SFE
66       COP = CFE
67       HPX = HPS*COP
68       HPY = HPS*SIP
69       IF (LT .EQ. 0) GO TO 2
70       WRITE(LUNOUT,3) HPS,HPX,HPY,HPZ,HE,HMA,LL
71     3 FORMAT(1H0,19HPS,PX,PY,PZ,E,M,LL=,6F8.4,I3)
72     2 CONTINUE
73       RETURN
74       END