source: trunk/MagicSoft/Simulation/Corsika/Mmcs/vapor.f@ 10083

Last change on this file since 10083 was 286, checked in by harald, 25 years ago
This is the start point for further developments of the Magic Monte Carlo Simulation written by Jose Carlos Gonzales. Now it is under control of one CVS repository for the whole collaboration. Everyone should use this CVS repository for further developments.
File size: 7.8 KB
Line 
1 SUBROUTINE VAPOR(MAPROJ,INEW,JFIN,ITYP,PFRX,PFRY)
2
3C-----------------------------------------------------------------------
4C (E)VAPOR(ATION OF NUCLEONS AND ALPHA PARTICLES FROM FRAGMENT)
5C
6C TREATES THE REMAINING UNFRAGMENTED NUCLEUS
7C EVAPORATION FOLLOWING CAMPI APPROXIMATION
8C SEE: X. CAMPI AND J. HUEFNER, PHYS.REV. C24 (1981) 2199
9C AND J.J. GAIMARD, THESE UNIVERSITE PARIS 7, (1990)
10C THIS SUBROUTINE IS CALLED FROM SDPM AND VSTORE
11C
12C ARGUMENTS INPUT:
13C MAPROJ = NUMBER OF NUCLEONS OF PROJECTILE
14C INEW = PARTICLE TYPE OF SPECTATOR FRAGMENT
15C ARGUMENTS OUTPUT:
16C JFIN = NUMBER OF FRAGMENTS
17C ITYP(1:JFIN) = NATURE (PARTICLE CODE) OF FRAGMENTS (GEANT)
18C PFRX(1:JFIN) = TRANSVERSE MOMENTUM OF FRAGMENTS IN X-DIRECTION
19C PFRY(1:JFIN) = TRANSVERSE MOMENTUM OF FRAGMENTS IN Y-DIRECTION
20C
21C DESIGN : D. HECK IK3 FZK KARLSRUHE
22C-----------------------------------------------------------------------
23
24 IMPLICIT NONE
25*KEEP,CONST.
26 COMMON /CONST/ PI,PI2,OB3,TB3,ENEPER
27 DOUBLE PRECISION PI,PI2,OB3,TB3,ENEPER
28*KEEP,DPMFLG.
29 COMMON /DPMFLG/ NFLAIN,NFLDIF,NFLPI0,NFLCHE,NFLPIF,NFRAGM
30 INTEGER NFLAIN,NFLDIF,NFLPI0,NFLCHE,NFLPIF,NFRAGM
31*KEEP,RANDPA.
32 COMMON /RANDPA/ FAC,U1,U2,RD,NSEQ,ISEED,KNOR
33 DOUBLE PRECISION FAC,U1,U2
34 REAL RD(3000)
35 INTEGER ISEED(103,10),NSEQ
36 LOGICAL KNOR
37*KEEP,RUNPAR.
38 COMMON /RUNPAR/ FIXHEI,THICK0,HILOECM,HILOELB,
39 * STEPFC,NRRUN,NSHOW,PATAPE,MONIIN,
40 * MONIOU,MDEBUG,NUCNUC,
41 * CETAPE,
42 * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
43 * N1STTR,MDBASE,
44 * DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
45 * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
46 * ,GHEISH,GHESIG
47 COMMON /RUNPAC/ DSN,HOST,USER
48 DOUBLE PRECISION FIXHEI,THICK0,HILOECM,HILOELB
49 REAL STEPFC
50 INTEGER NRRUN,NSHOW,PATAPE,MONIIN,MONIOU,MDEBUG,NUCNUC,
51 * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
52 * N1STTR,MDBASE
53 INTEGER CETAPE
54 CHARACTER*79 DSN
55 CHARACTER*20 HOST,USER
56
57 LOGICAL DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
58 * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
59 * ,GHEISH,GHESIG
60*KEND.
61
62 DOUBLE PRECISION PFR(60),PFRX(60),PFRY(60)
63 DOUBLE PRECISION AFIN,AGLH,APRF,BGLH,EEX,PHIFR,RANNOR,SPFRX,SPFRY
64 INTEGER ITYP(60),IARM,INEW,ITYPRM,INRM,IS,IZRM,JC,JFIN,
65 * K,L,LS,MAPROJ,MF,NFIN,NINTA,NNUC,NPRF,NSTEP
66 EXTERNAL RANNOR
67C-----------------------------------------------------------------------
68
69 IF (DEBUG) WRITE(MDEBUG,*) 'VAPOR : MAPROJ,INEW = ',MAPROJ,INEW
70
71 ITYPRM = INEW
72 NPRF = INEW/100
73 NINTA = MAPROJ - NPRF
74 IF ( NINTA .EQ. 0 ) THEN
75C NO NUCLEON HAS INTERACTED
76 JFIN = 1
77 PFR(1) = 0.D0
78 ITYP(1) = INEW
79 IF (DEBUG) WRITE(MDEBUG,*) 'VAPOR : JFIN,NINTA= ',JFIN,NINTA
80 RETURN
81 ENDIF
82
83C EXCITATION ENERGY EEX OF PREFRAGMENT
84C SEE: J.J. GAIMARD, THESE UNIVERSITE PARIS 7, (1990), CHPT. 4.2
85 EEX = 0.D0
86 CALL RMMAR(RD,2*NINTA,1)
87 DO 22 L = 1,NINTA
88 IF ( RD(NINTA+L) .LT. RD(L) ) RD(L) = 1. - RD(L)
89 EEX = EEX + RD(L)
90 22 CONTINUE
91C DEPTH OF WOODS-SAXON POTENTIAL TO FERMI SURFACE IS 0.040 GEV
92 IF (DEBUG) WRITE(MDEBUG,*)'VAPOR : EEX = ',SNGL(EEX*0.04D0),' GEV'
93C EVAPORATION: EACH EVAPORATION STEP NEEDS ABOUT 0.020 GEV, THEREFORE
94C NSTEP IS EEX * 0.04/0.02 = EEX * 2.
95 NSTEP = INT(EEX*2.D0)
96
97 IF ( NSTEP .LE. 0 ) THEN
98C EXCITATION ENERGY TOO SMALL, NO EVAPORATION
99 JFIN = 1
100 PFR(1) = 0.D0
101 ITYP(1) = INEW
102 IF (DEBUG) WRITE(MDEBUG,*) 'VAPOR : JFIN,EEX = ',JFIN,SNGL(EEX)
103 RETURN
104 ENDIF
105
106C AFIN IS ATOMIC NUMBER OF FINAL NUCLEUS
107 APRF = FLOAT(NPRF)
108 AFIN = APRF - 1.6D0 * FLOAT(NSTEP)
109 NFIN = MAX( INT(AFIN+0.5D0), 0 )
110C CORRESPONDS TO DEFINITION; FRAGMENTATION-EVAPORATION
111C CONVOLUTION EMU07 /MODEL ABRASION EVAPORATION (JNC FZK APRIL 94)
112C NNUC IS NUMBER OF EVAPORATING NUCLEONS
113 NNUC = NPRF - NFIN
114 IF (DEBUG) WRITE(MDEBUG,*) 'VAPOR : NFIN,NNUC = ',NFIN,NNUC
115 JC = 0
116
117 IF ( NNUC .LE. 0 ) THEN
118C NO EVAPORATION
119 JFIN = 1
120 PFR(1) = 0.D0
121 ITYP(1) = INEW
122 RETURN
123
124 ELSEIF ( NNUC .GE. 4 ) THEN
125C EVAPORATION WITH FORMATION OF ALPHA PARTICLES POSSIBLE
126C IARM, IZRM, INRM ARE NUMBER OF NUCLEONS, PROTONS, NEUTRONS OF
127C REMAINDER
128 DO 31 LS = 1,NSTEP
129 IARM = ITYPRM/100
130 IF ( IARM .LE. 0 ) GOTO 100
131 IZRM = MOD(ITYPRM,100)
132 INRM = IARM - IZRM
133 JC = JC + 1
134 CALL RMMAR(RD,2,1)
135 IF ( RD(1).LT.0.2 .AND. IZRM.GE.2 .AND. INRM.GE.2 ) THEN
136 ITYP(JC) = 402
137 NNUC = NNUC - 4
138 ITYPRM = ITYPRM - 402
139 ELSE
140 IF ( RD(2)*(IZRM+INRM) .LT. IZRM ) THEN
141 ITYP(JC) = 14
142 ITYPRM = ITYPRM - 101
143 ELSE
144 ITYP(JC) = 13
145 ITYPRM = ITYPRM - 100
146 ENDIF
147 NNUC = NNUC - 1
148 ENDIF
149 IF ( NNUC .LE. 0 ) GOTO 50
150 31 CONTINUE
151 ENDIF
152
153 IF ( NNUC .LT. 4 ) THEN
154C EVAPORATION WITHOUT FORMATION OF ALPHA PARTICLES
155 CALL RMMAR(RD,NNUC,1)
156 DO 32 IS = 1,NNUC
157 IARM = ITYPRM/100
158 IF ( IARM .LE. 0 ) GOTO 100
159 IZRM = MOD(ITYPRM,100)
160 JC = JC + 1
161 IF ( RD(IS)*IARM .LT. IZRM ) THEN
162 ITYP(JC) = 14
163 ITYPRM = ITYPRM - 101
164 ELSE
165 ITYP(JC) = 13
166 ITYPRM = ITYPRM - 100
167 ENDIF
168 32 CONTINUE
169 ENDIF
170
171 50 CONTINUE
172 JC = JC + 1
173 IF ( ITYPRM .GT. 101 ) THEN
174 ITYP(JC) = ITYPRM
175 ELSEIF ( ITYPRM .EQ. 101 ) THEN
176 ITYP(JC) = 14
177 ELSEIF ( ITYPRM .EQ. 100 ) THEN
178 ITYP(JC) = 13
179 ELSE
180 JC = JC - 1
181 IF ( ITYPRM .NE. 0 ) WRITE(MONIOU,*)
182 * 'VAPOR : ILLEGAL PARTICLE ITYPRM =',ITYPRM
183 ENDIF
184
185 100 JFIN = JC
186 IF (DEBUG) WRITE(MDEBUG,*) 'VAPOR : NO ITYP PFR'
187 IF ( NFRAGM .EQ. 2 ) THEN
188C EVAPORATION WITH PT AFTER PARAMETRIZED JACEE DATA
189 DO 150 MF = 1,JFIN
190 PFR(MF) = RANNOR(0.088D0,0.044D0)
191 IF (DEBUG) WRITE(MDEBUG,*) MF,ITYP(MF),SNGL(PFR(MF))
192 150 CONTINUE
193 ELSEIF ( NFRAGM .EQ. 3 ) THEN
194C EVAPORATION WITH PT AFTER GOLDHABER'S MODEL (PHYS.LETT.53B(1974)306)
195 DO 160 MF = 1,JFIN
196 K = MAX( 1, ITYP(MF)/100 )
197 BGLH = K * (MAPROJ - K) / FLOAT(MAPROJ-1)
198C THE VALUE 0.103 [GEV] IS SIGMA(0)=P(FERMI)/SQRT(5.)
199* AGLH = 0.103D0 * SQRT( BGLH )
200C THE VALUE 0.090 [GEV] IS EXPERIMENTALLY DETERMINED SIGMA(0)
201 AGLH = 0.090D0 * SQRT( BGLH )
202 PFR(MF) = RANNOR(0.D0,AGLH)
203 IF (DEBUG) WRITE(MDEBUG,*) MF,ITYP(MF),SNGL(PFR(MF))
204 160 CONTINUE
205 ELSE
206C EVAPORATION WITHOUT TRANSVERSE MOMENTUM
207 DO 165 MF = 1,JFIN
208 PFR(MF) = 0.D0
209 IF (DEBUG) WRITE(MDEBUG,*) MF,ITYP(MF),SNGL(PFR(MF))
210 165 CONTINUE
211 ENDIF
212C CALCULATE RESIDUAL TRANSVERSE MOMENTUM
213 SPFRX = 0.D0
214 SPFRY = 0.D0
215 CALL RMMAR(RD,JFIN,1)
216 DO 170 MF = 1,JFIN
217 PHIFR = PI * RD(MF)
218 PFRX(MF) = PFR(MF) * COS(PHIFR)
219 PFRY(MF) = PFR(MF) * SIN(PHIFR)
220 SPFRY = SPFRY + PFRY(MF)
221 SPFRX = SPFRX + PFRX(MF)
222 170 CONTINUE
223C CORRECT ALL TRANSVERSE MOMENTA FOR MOMENTUM CONSERVATION
224 SPFRX = SPFRX / JFIN
225 SPFRY = SPFRY / JFIN
226 DO 180 MF = 1,JFIN
227 PFRX(MF) = PFRX(MF) - SPFRX
228 PFRY(MF) = PFRY(MF) - SPFRY
229 180 CONTINUE
230
231 IF (DEBUG) WRITE(MDEBUG,*) 'VAPOR : NINTA,JFIN= ',NINTA,JFIN
232 RETURN
233 END
Note: See TracBrowser for help on using the repository browser.