This commit was generated by cvs2svn to compensate for changes in r2,
[u/mrichter/AliRoot.git] / GEANT321 / peanut / rstsel.F
1 *
2 * $Id$
3 *
4 * $Log$
5 * Revision 1.1.1.1  1995/10/24 10:22:03  cernlib
6 * Geant
7 *
8 *
9 #include "geant321/pilot.h"
10 *CMZ :  3.21/02 29/03/94  15.41.46  by  S.Giani
11 *-- Author :
12 *$ CREATE RSTSEL.FOR
13 *COPY RSTSEL
14 *
15 *=== rstsel ===========================================================*
16 *
17       SUBROUTINE RSTSEL ( JPROJ )
18  
19 #include "geant321/dblprc.inc"
20 #include "geant321/dimpar.inc"
21 #include "geant321/iounit.inc"
22 *
23 *----------------------------------------------------------------------*
24 *----------------------------------------------------------------------*
25 *
26 #include "geant321/balanc.inc"
27 #include "geant321/nucdat.inc"
28 #include "geant321/nucgeo.inc"
29 #include "geant321/part.inc"
30 #include "geant321/parevt.inc"
31 #include "geant321/resnuc.inc"
32 *
33       REAL RNDM(1)
34       SAVE ABTAR , ZZTAR , PACORE, PASKIN,
35      &     PAHALO, PACHLP, XHLP  , AUSFL , ZUSFL , ONUSFL, KPROJ ,
36      &     ITFRMI, ITFRM2
37 *
38       KPROJ  = JPROJ
39       ABTAR  = IBTAR
40       ZZTAR  = ICHTAR
41       AUSFL  = ABTAR
42       ZUSFL  = ZZTAR
43       ONUSFL = 0.D+00
44       CALL RACO ( CXIMPC, CYIMPC, CZIMPC )
45       UBIMPC = CXIMPC
46       VBIMPC = CYIMPC
47       WBIMPC = CZIMPC
48       PACORE = FRHINC (RADIU0) / ABTAR
49       PASKIN = FRHINC (RADIU1) / ABTAR - PACORE
50       PAHALO = 1.D+00 - PACORE - PASKIN
51       GO TO 500
52       ENTRY RSTNXT
53   500 CONTINUE
54       LABRST = .TRUE.
55       IF ( KPROJ .EQ. 2 ) THEN
56          NTARGT = 1
57          STOP 'pbar_rest'
58       ELSE IF ( KPROJ .EQ. 9 ) THEN
59          NTARGT = 1
60          STOP 'nbar_rest'
61       ELSE IF ( KPROJ .EQ. 14 ) THEN
62          CALL GRNDM(RNDM,1)
63          IF ( RNDM (1) .LT. RADPIM * ( ABTAR + RDPMHL ) / ABTAR )
64      &      THEN
65             NTARGT = 1
66             KNUCIM = 1
67             ITFRMI = 1
68             IPRTYP = - ( 14 * 10 + KNUCIM )
69          ELSE
70             NTARGT = 2
71             KNUCIM = 1
72             ITFRMI = 1
73             PRCOL0 = APMRST
74             PRCOLP = PRCOL0 * ( ZUSFL - ONUSFL )
75             PRCOLP = PRCOLP / ( PRCOLP + ( 1.D+00 - PRCOL0 ) * ( AUSFL
76      &             - ZUSFL ) )
77             CALL GRNDM(RNDM,1)
78             IF ( RNDM (1) .LT. PRCOLP ) THEN
79                KNUCI2 = 1
80                ITFRM2 = 1
81             ELSE
82                KNUCI2 = 8
83                ITFRM2 = 2
84             END IF
85             IPRTYP = - ( 14 * 100 + KNUCIM * 10 + KNUCI2 )
86          END IF
87       ELSE IF ( KPROJ .EQ. 15 ) THEN
88          NTARGT = 1
89          NTARGT = 2
90          STOP 'K-_rest'
91       ELSE
92          STOP '???_rest'
93       END IF
94       CALL GRNDM(RNDM,1)
95       RNDMAB = RNDM (1)
96       IF ( RNDMAB .LT. PACORE ) THEN
97          RNDMAB = RNDMAB / PACORE
98          RADSAM = RADIU0
99          RHOSAM = RHOCEN
100          XHLP   = 0.D+00
101       ELSE IF ( RNDMAB - PACORE .LT. PASKIN ) THEN
102          RNDMAB = ( RNDMAB - PACORE ) / PASKIN
103          RADSAM = RADIU1
104          RHOSAM = RHOCOR
105          XHLP   = RADIU0 / RADIU1
106       ELSE
107          RNDMAB = ( RNDMAB - PACORE - PASKIN ) / PAHALO
108          RHOSAM = RHOSKN
109          RADSAM = RADTOT
110          XHLP   = RADIU1 / RADTOT
111       END IF
112       XHLP3  = XHLP * XHLP * XHLP
113  1000 CONTINUE
114          BIMPCT = RADSAM * ( RNDMAB * ( 1.D+00 - XHLP3 ) + XHLP3 )
115      &          **0.3333333333333333D+00
116          RHOIMP = FRHONC (BIMPCT)
117          CALL GRNDM(RNDM,1)
118          RNDMAB = RNDM (1)
119          RHORAT = RHOIMP / RHOSAM
120          IF ( RNDMAB .GE. RHORAT ) THEN
121             RNDMAB = ( RNDMAB - RHORAT ) / ( 1.D+00 - RHORAT )
122             GO TO 1000
123          END IF
124       PFRIMP = FPFRNC ( RHOIMP, ITFRMI )
125       EKFIMP = FEKFNC ( PFRIMP, ITFRMI )
126       RHOIMT = RHOIMP
127       IF ( NTARGT .EQ. 2 ) THEN
128          IF ( RHOIMP .GE. RHOCEN ) THEN
129             PFRIM2 = PFRCEN (ITFRM2)
130             EKFIM2 = EKFCEN (ITFRM2)
131             EKFPRO = EKFCEN (IPWELL)
132             PFRPRO = PFRCEN (IPWELL)
133          ELSE IF ( ITFRM2 .NE. ITFRMI ) THEN
134             PFRIM2 = PFRIMP * PFRCEN (ITFRM2) / PFRCEN (ITFRMI)
135             EKFIM2 = SQRT ( PFRIM2 * PFRIM2 + AMNUSQ (ITFRM2) )
136      &             - AMNUCL (ITFRM2)
137             IF ( IPWELL .EQ. ITFRMI ) THEN
138                EKFPRO = EKFIMP
139                PFRPRO = PFRIMP
140             ELSE
141                PFRPRO = PFRIM2
142                EKFPRO = EKFIM2
143             END IF
144          ELSE
145             PFRIM2 = PFRIMP
146             EKFIM2 = EKFIMP
147             IF ( IPWELL .EQ. ITFRMI ) THEN
148                PFRPRO = PFRIMP
149                EKFPRO = EKFIMP
150             ELSE
151                PFRPRO = PFRIMP * PFRCEN (IPWELL) / PFRCEN (ITFRMI)
152                EKFPRO = SQRT ( PFRPRO * PFRPRO + AMNUSQ (IPWELL) )
153      &                - AMNUCL (IPWELL)
154             END IF
155          END IF
156       ELSE
157          IF ( RHOIMP .GE. RHOCEN ) THEN
158             PFRPRO = PFRCEN (IPWELL)
159             EKFPRO = EKFCEN (IPWELL)
160          ELSE IF ( IPWELL .EQ. ITFRMI ) THEN
161             PFRPRO = PFRIMP
162             EKFPRO = EKFIMP
163          ELSE
164             PFRPRO = PFRIMP * PFRCEN (IPWELL) / PFRCEN (ITFRMI)
165             EKFPRO = SQRT ( PFRPRO * PFRPRO + AMNUSQ (IPWELL) )
166      &             - AMNUCL (IPWELL)
167          END IF
168       END IF
169       VPRWLL = WLLRED * ( EKFPRO + BNDNUC )
170       BIMPTR = BIMPCT
171       BIMMEM = BIMPTR
172       RIMPCT = BIMPCT
173       RIMPTR = BIMPTR
174       XBIMPC = UBIMPC * BIMPCT
175       YBIMPC = VBIMPC * BIMPCT
176       ZBIMPC = WBIMPC * BIMPCT
177       XIMPCT = XBIMPC
178       YIMPCT = YBIMPC
179       ZIMPCT = ZBIMPC
180       XIMPTR = XIMPCT
181       YIMPTR = YIMPCT
182       ZIMPTR = ZIMPCT
183       EKEWLL = EKECON + VPRWLL
184       PPRWLL = SQRT ( EKEWLL * ( EKEWLL + 2.D+00 * AM (KPROJ) ) )
185       PXPROJ = PPRWLL * CXIMPC
186       PYPROJ = PPRWLL * CYIMPC
187       PZPROJ = PPRWLL * CZIMPC
188       PNFRMI = PFNCLV ( ITFRMI, .TRUE. )
189       IF ( PNFRMI .LT. -100.D+00 ) GO TO 500
190       CALL RACO ( PXFERM, PYFERM, PZFERM )
191       PXFERM = PXFERM * PNFRMI
192       PYFERM = PYFERM * PNFRMI
193       PZFERM = PZFERM * PNFRMI
194       EKFIMP = MAX ( EKFERM, EKFIMP )
195       IF ( NTARGT .EQ. 2 ) THEN
196          PNFRM2 = PFNCLV ( ITFRM2, .FALSE. )
197          IF ( PNFRM2 .LT. -100.D+00 ) GO TO 500
198          CALL RACO ( PXFER2, PYFER2, PZFER2 )
199          PXFER2 = PXFER2 * PNFRM2
200          PYFER2 = PYFER2 * PNFRM2
201          PZFER2 = PZFER2 * PNFRM2
202          EKFIM2 = MAX ( EKFER2, EKFIM2 )
203          ERES   = EKEWLL + AM (KPROJ) + AM (KNUCIM) + EKFERM
204      &          + AM (KNUCI2) + EKFER2
205          PXRES  = PXPROJ + PXFERM + PXFER2
206          PYRES  = PYPROJ + PYFERM + PYFER2
207          PZRES  = PZPROJ + PZFERM + PZFER2
208          PTRES2 = PXRES**2 + PYRES**2 + PZRES**2
209       ELSE
210          ERES   = EKEWLL + AM (KPROJ) + AM (KNUCIM) + EKFERM
211          PXRES  = PXPROJ + PXFERM
212          PYRES  = PYPROJ + PYFERM
213          PZRES  = PZPROJ + PZFERM
214          PTRES2 = PXRES**2 + PYRES**2 + PZRES**2
215       END IF
216       RETURN
217 *=== End of subroutine rstsel =========================================*
218       END