source: trunk/MagicSoft/Simulation/Corsika/Mmcs/sdpm.f@ 6495

Last change on this file since 6495 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: 19.7 KB
Line 
1 SUBROUTINE SDPM
2
3C-----------------------------------------------------------------------
4C S(TARTING) D(UAL) P(ARTON) M(ODEL)
5C
6C THIS ROUTINE DETERMINES THE TARGET NUCLEUS.
7C IT CALLS ALSO THE VARIOUS INTERACTION MODELS.
8C FOR HDPM, THIS ROUTINE LOOKS, HOW MANY NUCLEONS INTERACT AND WHICH
9C RESIDUAL FRAGMENT OF THE PROJECTILE NUCLEUS REMAINS.
10C THIS SUBROUTINE IS CALLED FROM NUCINT AND PIGEN
11C
12C REDESIGN: D. HECK IK3 FZK KARLSRUHE
13C-----------------------------------------------------------------------
14
15 IMPLICIT NONE
16*KEEP,AIR.
17 COMMON /AIR/ COMPOS,PROBTA,AVERAW,AVOGAD
18 DOUBLE PRECISION COMPOS(3),PROBTA(3),AVERAW,AVOGAD
19*KEEP,CONST.
20 COMMON /CONST/ PI,PI2,OB3,TB3,ENEPER
21 DOUBLE PRECISION PI,PI2,OB3,TB3,ENEPER
22*KEEP,DPMFLG.
23 COMMON /DPMFLG/ NFLAIN,NFLDIF,NFLPI0,NFLCHE,NFLPIF,NFRAGM
24 INTEGER NFLAIN,NFLDIF,NFLPI0,NFLCHE,NFLPIF,NFRAGM
25*KEEP,GENER.
26 COMMON /GENER/ GEN,ALEVEL
27 DOUBLE PRECISION GEN,ALEVEL
28*KEEP,INTER.
29 COMMON /INTER/ AVCH,AVCH3,DC0,DLOG,DMLOG,ECMDIF,ECMDPM,ELAB,
30 * FNEUT,FNEUT2,GNU,PLAB,POSC2,POSC3,POSN2,POSN3,
31 * RC3TO2,S,SEUGF,SEUGP,SLOG,SLOGSQ,SMLOG,
32 * WIDC2,WIDC3,WIDN2,WIDN3,YCM,YY0,ZN,
33 * IDIF,ITAR
34 DOUBLE PRECISION AVCH,AVCH3,DC0,DLOG,DMLOG,ECMDIF,ECMDPM,ELAB,
35 * FNEUT,FNEUT2,GNU,PLAB,POSC2,POSC3,POSN2,POSN3,
36 * RC3TO2,S,SEUGF,SEUGP,SLOG,SLOGSQ,SMLOG,
37 * WIDC2,WIDC3,WIDN2,WIDN3,YCM,YY0,ZN
38 INTEGER IDIF,ITAR
39*KEEP,ISTA.
40 COMMON /ISTA/ IFINET,IFINNU,IFINKA,IFINPI,IFINHY
41 INTEGER IFINET,IFINNU,IFINKA,IFINPI,IFINHY
42*KEEP,MULT.
43 COMMON /MULT/ EKINL,MSMM,MULTMA,MULTOT
44 DOUBLE PRECISION EKINL
45 INTEGER MSMM,MULTMA(37,13),MULTOT(37,13)
46*KEEP,NCSNCS.
47 COMMON /NCSNCS/ SIGN30,SIGN45,SIGN60,SIGO30,SIGO45,SIGO60,
48 * SIGA30,SIGA45,SIGA60,PNOA30,PNOA45,PNOA60,
49 * SIG30A,SIG45A,SIG60A
50 DOUBLE PRECISION SIGN30(56),SIGN45(56),SIGN60(56),
51 * SIGO30(56),SIGO45(56),SIGO60(56),
52 * SIGA30(56),SIGA45(56),SIGA60(56),
53 * PNOA30(1540,3),PNOA45(1540,3),PNOA60(1540,3),
54 * SIG30A(56),SIG45A(56),SIG60A(56)
55*KEEP,PAM.
56 COMMON /PAM/ PAMA,SIGNUM
57 DOUBLE PRECISION PAMA(6000),SIGNUM(6000)
58*KEEP,PARPAR.
59 COMMON /PARPAR/ CURPAR,SECPAR,PRMPAR,OUTPAR,C,
60 * E00,E00PN,PTOT0,PTOT0N,THICKH,ITYPE,LEVL
61 DOUBLE PRECISION CURPAR(14),SECPAR(14),PRMPAR(14),OUTPAR(14),
62 * C(50),E00,E00PN,PTOT0,PTOT0N,THICKH
63 INTEGER ITYPE,LEVL
64*KEEP,PARPAE.
65 DOUBLE PRECISION GAMMA,COSTHE,PHI,H,T,X,Y,CHI,BETA,GCM,ECM
66 EQUIVALENCE (CURPAR(2),GAMMA), (CURPAR(3),COSTHE),
67 * (CURPAR(4), PHI ), (CURPAR(5), H ),
68 * (CURPAR(6), T ), (CURPAR(7), X ),
69 * (CURPAR(8), Y ), (CURPAR(9), CHI ),
70 * (CURPAR(10),BETA), (CURPAR(11),GCM ),
71 * (CURPAR(12),ECM )
72*KEEP,RANDPA.
73 COMMON /RANDPA/ FAC,U1,U2,RD,NSEQ,ISEED,KNOR
74 DOUBLE PRECISION FAC,U1,U2
75 REAL RD(3000)
76 INTEGER ISEED(103,10),NSEQ
77 LOGICAL KNOR
78*KEEP,RANGE.
79 COMMON /RANGE/ CC
80 DOUBLE PRECISION CC(20)
81*KEEP,REST.
82 COMMON /REST/ CONTNE,TAR,LT
83 DOUBLE PRECISION CONTNE(3),TAR
84 INTEGER LT
85*KEEP,RUNPAR.
86 COMMON /RUNPAR/ FIXHEI,THICK0,HILOECM,HILOELB,
87 * STEPFC,NRRUN,NSHOW,PATAPE,MONIIN,
88 * MONIOU,MDEBUG,NUCNUC,
89 * CETAPE,
90 * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
91 * N1STTR,MDBASE,
92 * DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
93 * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
94 * ,GHEISH,GHESIG
95 COMMON /RUNPAC/ DSN,HOST,USER
96 DOUBLE PRECISION FIXHEI,THICK0,HILOECM,HILOELB
97 REAL STEPFC
98 INTEGER NRRUN,NSHOW,PATAPE,MONIIN,MONIOU,MDEBUG,NUCNUC,
99 * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
100 * N1STTR,MDBASE
101 INTEGER CETAPE
102 CHARACTER*79 DSN
103 CHARACTER*20 HOST,USER
104
105 LOGICAL DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
106 * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
107 * ,GHEISH,GHESIG
108*KEEP,SIGM.
109 COMMON /SIGM/ SIGMA,SIGANN,SIGAIR,FRACTN,FRCTNO
110 DOUBLE PRECISION SIGMA,SIGANN,SIGAIR,FRACTN,FRCTNO
111*KEEP,VKIN.
112 COMMON /VKIN/ BETACM
113 DOUBLE PRECISION BETACM
114*KEEP,VENUS.
115 COMMON /VENUS/ ISH0,IVERVN,MTAR99,FVENUS,FVENSG
116 INTEGER ISH0,IVERVN,MTAR99
117 LOGICAL FVENUS,FVENSG
118*KEND.
119
120 DOUBLE PRECISION PFRX(60),PFRY(60)
121 DOUBLE PRECISION COSTET,EA,P,PHIV,PTM,PT2,
122 * SIGMAA,SIGMAN,SIGMAO,SIG45,S45SQ,S4530
123 DOUBLE PRECISION CGHSIG,EKIN
124 EXTERNAL CGHSIG
125 INTEGER ITYP(60),I,IA,IANEW,INACTA,INACTZ,INDEX,INEUTR,
126 * IZ,IZNEW,J,JFIN,KNEW,L,LL,NPRPRO,NNEPRO
127C-----------------------------------------------------------------------
128
129 IF ( DEBUG ) WRITE(MDEBUG,444) (CURPAR(I),I=1,9)
130 444 FORMAT(' SDPM : CURPAR=',1P,9E10.3)
131
132C IA IS MASS NUMBER OF PROJECTILE
133 IA = ITYPE / 100
134 IF ( IA .GT. 56 ) THEN
135 WRITE(MONIOU,*) 'SDPM : NOT FORESEEN PARTICLE TYPE=',ITYPE
136 STOP
137 ENDIF
138
139C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
140C TREATMENT OF GAMMAS COMING FROM EGS4 (PIGEN)
141 IF ( ITYPE .EQ. 1 ) THEN
142C RATIOS OF CROSS SECTIONS GO LIKE A**0.91
143 FRACTN = COMPOS(1) * 11.04019D0
144 FRCTNO = FRACTN + COMPOS(2) * 12.46663D0
145 SIGAIR = FRCTNO + COMPOS(3) * 28.69952D0
146C TARGET IS CHOSEN AT RANDOM
147 CALL RMMAR( RD,1,1 )
148 IF ( RD(1)*SIGAIR .LE. FRACTN ) THEN
149C INTERACTION WITH NITROGEN
150 LT = 1
151 TAR = 14.D0
152 ELSEIF ( RD(1)*SIGAIR .LE. FRCTNO ) THEN
153C INTERACTION WITH OXYGEN
154 LT = 2
155 TAR = 16.D0
156 ELSE
157C INTERACTION WITH ARGON
158 LT = 3
159 TAR = 40.D0
160 ENDIF
161
162C GAMMAS ARE TREATED BY VENUS, IF SUFFICIENT ENERGY
163 IF ( FVENUS .AND. CURPAR(2) .GT. HILOELB ) THEN
164 CALL VENLNK
165 ELSE
166 CALL HDPM
167 ENDIF
168
169C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
170C NORMAL HADRON PROJECTILE
171 ELSEIF ( ITYPE .LT. 100 ) THEN
172
173C WITH WHAT KIND OF TARGET DOES PROJECTILE INTERACT?
174 IF ( FIXTAR ) THEN
175C TARGET OF FIRST INTERACTION IS FIXED
176 LT = N1STTR
177 IF ( N1STTR .EQ. 1 ) THEN
178 TAR = 14.D0
179 ELSEIF ( N1STTR .EQ. 2 ) THEN
180 TAR = 16.D0
181 ELSE
182 TAR = 40.D0
183 ENDIF
184 FIXTAR = .FALSE.
185 ELSE
186C TARGET IS CHOSEN AT RANDOM ACCORDING TO CROSS SECTION
187C SIGMA IS ENERGY DEPENDENT INELASTIC NUCLEON-NUCLEON CROSS SECTION
188C AND IS SET IN BOX2
189C AUXIL. QUANTITIES FOR INTERPOLATION
190 SIG45 = SIGMA - 45.D0
191 S45SQ = SIG45**2 / 450.D0
192 S4530 = SIG45 / 30.D0
193C INELASTIC CROSS SECTIONS FOR PROJECTICLE WITH MASS NUMBER 1
194 SIGMAN = (1.D0 - 2.D0 * S45SQ) * SIGN45(1)
195 * +(S45SQ - S4530) * SIGN30(1)
196 * +(S45SQ + S4530) * SIGN60(1)
197 FRACTN = COMPOS(1) * SIGMAN
198 SIGMAO = (1.D0 - 2.D0 * S45SQ) * SIGO45(1)
199 * +(S45SQ - S4530) * SIGO30(1)
200 * +(S45SQ + S4530) * SIGO60(1)
201 FRCTNO = FRACTN + COMPOS(2) * SIGMAO
202 SIGMAA = (1.D0 - 2.D0 * S45SQ) * SIGA45(1)
203 * +(S45SQ - S4530) * SIGA30(1)
204 * +(S45SQ + S4530) * SIGA60(1)
205C INELASTIC CROSS SECTIONS OF AIR FOR PROJECTILE WITH MASS NUMBER 1
206 SIGAIR = FRCTNO + COMPOS(3)*SIGMAA
207 333 CONTINUE
208 CALL RMMAR( RD,1,1 )
209 IF(DEBUG)WRITE(MDEBUG,*)'SDPM : FRACTN=',SNGL(FRACTN),
210 * 'FRCTNO=',SNGL(FRCTNO),'RD=',RD(1)
211 IF ( RD(1)*SIGAIR .LE. FRACTN ) THEN
212C INTERACTION WITH NITROGEN
213 LT = 1
214 TAR = 14.D0
215 ELSEIF ( RD(1)*SIGAIR .LE. FRCTNO ) THEN
216C INTERACTION WITH OXYGEN
217 LT = 2
218 TAR = 16.D0
219 ELSE
220C INTERACTION WITH ARGON
221 LT = 3
222 TAR = 40.D0
223 ENDIF
224 ENDIF
225
226 IF ( FVENUS ) THEN
227C MESONS, NUCLEONS AND STRANGE BARYONS ARE TREATED BY VENUS (JAN 95)
228 IF ( (ITYPE .GE. 7 .AND. ITYPE .LE. 16) .OR.
229 * (ITYPE .GE. 18 .AND. ITYPE .LE. 32) )THEN
230 CALL VENLNK
231 ELSE
232 CALL HDPM
233 ENDIF
234 ELSE
235 CALL HDPM
236 ENDIF
237
238C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
239C HEAVY PRIMARY INCIDENT WITH IA NUCLEONS
240 ELSEIF ( IA .LE. 56 ) THEN
241
242 IZ = MOD(ITYPE,100)
243C WITH WHAT KIND OF TARGET DOES PROJECTILE INTERACT?
244 IF ( FIXTAR ) THEN
245C TARGET OF FIRST INTERACTION IS FIXED
246 LT = N1STTR
247 IF ( N1STTR .EQ. 1 ) THEN
248 TAR = 14.D0
249 ELSEIF ( N1STTR .EQ. 2 ) THEN
250 TAR = 16.D0
251 ELSE
252 TAR = 40.D0
253 ENDIF
254 FIXTAR = .FALSE.
255 CALL RMMAR( RD,2,1 )
256 ELSE
257C ONLY INELASTIC INTERACTIONS WITH HEAVY PRIMARY/FRAGMENT
258C SIGMA IS ENERGY DEPENDENT INELASTIC NUCLEON-NUCLEON CROSS SECTION
259C AND IS SET IN BOX2
260C AUXIL. QUANTITIES FOR INTERPOLATION
261 SIG45 = SIGMA - 45.D0
262 S45SQ = SIG45**2 / 450.D0
263 S4530 = SIG45 / 30.D0
264C INELASTIC CROSS SECTIONS FOR PROJECTICLE WITH MASS NUMBER IA
265 SIGMAN = (1.D0 - 2.D0 * S45SQ) * SIGN45(IA)
266 * +(S45SQ - S4530) * SIGN30(IA)
267 * +(S45SQ + S4530) * SIGN60(IA)
268 FRACTN = COMPOS(1) * SIGMAN
269 SIGMAO = (1.D0 - 2.D0 * S45SQ) * SIGO45(IA)
270 * +(S45SQ - S4530) * SIGO30(IA)
271 * +(S45SQ + S4530) * SIGO60(IA)
272 FRCTNO = FRACTN + COMPOS(2) * SIGMAO
273 SIGMAA = (1.D0 - 2.D0 * S45SQ) * SIGA45(IA)
274 * +(S45SQ - S4530) * SIGA30(IA)
275 * +(S45SQ + S4530) * SIGA60(IA)
276C INELASTIC CROSS SECTIONS OF AIR FOR PROJECTILE WITH MASS NUMBER IA
277 SIGAIR = FRCTNO +COMPOS(3)*SIGMAA
278C TARGET IS CHOSEN AT RANDOM
279 CALL RMMAR( RD,2,1 )
280 IF(DEBUG)WRITE(MDEBUG,*)'SDPM : FRACTN=',SNGL(FRACTN),
281 * 'FRCTNO=',SNGL(FRCTNO),'RD=',RD(1)
282 IF ( RD(1)*SIGAIR .LE. FRACTN ) THEN
283C INTERACTION WITH NITROGEN
284 LT = 1
285 TAR = 14.D0
286 ELSEIF ( RD(1)*SIGAIR .LE. FRCTNO ) THEN
287C INTERACTION WITH OXYGEN
288 LT = 2
289 TAR = 16.D0
290 ELSE
291C INTERACTION WITH ARGON
292 LT = 3
293 TAR = 40.D0
294 ENDIF
295 ENDIF
296C TREAT NUCLEUS BY VENUS, IF SELECTED AND ENERGY/NUCLEON HIGH ENOUGH
297 IF ( FVENUS .AND. PAMA(ITYPE)*GAMMA .GT. HILOELB*IA ) THEN
298 CALL VENLNK
299 RETURN
300 ENDIF
301
302C TREATMENT OF NUCLEUS-NUCLEUS INTERACTION IN HDPM BY SUPERPOSITION
303C
304C INDEX CALCULATION 1<I=<56 NUCLEONS IN PROJECTILE
305C 1<J<I INTERACTING NUCLEONS
306C P(I,I)=1 CUMULATIVE PROBABILITIES
307C P(I,J) ---> P( I*(I-3)*0.5+J+1 )
308C IZ IS NUMBER OF PROTONS IN PROJECTILE
309C LT IS INDEX FOR TARGET 1 = N, 2 = O, 3 = AR
310C INACTA IS NUMBER OF INTERACTING NUCLEONS
311C INACTZ IS NUMBER OF INTERACTING PROTONS
312
313C LOOK, HOW MANY NUCLEONS INTERACT
314 DO 100 J = 1,IA-1
315 INACTA = J
316 INDEX = IA * (IA-3) * 0.5 + 1 + J
317 P = ( 1.D0 - S45SQ *2.D0 ) * PNOA45(INDEX,LT)
318 * +( S45SQ - S4530 ) * PNOA30(INDEX,LT)
319 * +( S45SQ + S4530 ) * PNOA60(INDEX,LT)
320 IF ( RD(2) .LT. P ) GO TO 110
321 100 CONTINUE
322C ALL NUCLEONS INTERACT (INACTA EQUAL IA)
323 INACTA = INACTA + 1
324
325 110 CONTINUE
326 IANEW = IA - INACTA
327
328C REMAINING PROJECTILE WITH IANEW NUCLEONS
329 DO 120 L = 2,8
330 SECPAR(L) = CURPAR(L)
331 120 CONTINUE
332
333
334C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
335C PROJECTILE NUCLEUS FRAGMENTS COMPLETELY, WRITE SPECTATOR NUCLEONS
336C ONTO STACK
337 IF ( NFRAGM .EQ. 0 ) THEN
338C LOOK, HOW MANY PROTONS AND NEUTRONS ARE FORMED
339 IZNEW = IANEW / 2.15D0 + 0.7D0
340 INEUTR = IANEW - IZNEW
341 INACTZ = MAX( IZ-IZNEW, 0 )
342
343 IF ( IZNEW .GT. 0 ) THEN
344C PROTONS
345 SECPAR(1) = 14.D0
346 DO 300 L = 1,IZNEW
347 CALL TSTACK
348 300 CONTINUE
349 ENDIF
350 IF ( INEUTR .GT. 0 ) THEN
351C NEUTRONS
352 SECPAR(1) = 13.D0
353 DO 310 L = 1,INEUTR
354 CALL TSTACK
355 310 CONTINUE
356 ENDIF
357
358C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
359C NO FRAGMENTATION, BUT SUCCESSIVE ABRASION OF PROJECTILE NUCLEUS
360 ELSE
361 IF ( DEBUG ) WRITE( MDEBUG,111 ) TAR,INACTA,IANEW
362 111 FORMAT(' SDPM : TARGET=',F4.0,' INACTA=',I4,' IANEW=',I4)
363
364C ALL NUCLEONS INTERACT, NO RESIDUAL NUCLEUS
365 IF ( IANEW .EQ. 0 ) THEN
366 INACTZ = IZ
367 IF ( DEBUG ) WRITE(MDEBUG,554) (CURPAR(I),I=1,9)
368 554 FORMAT (' SDPM : CURPAR=',1P,9E10.3)
369 KNEW = 0
370
371C REMAINING NUCLEUS IS A NUCLEON
372 ELSEIF ( IANEW .EQ. 1 ) THEN
373 CALL RMMAR( RD,1,1 )
374 IZNEW = NINT(RD(1))
375 INACTZ = IZ - IZNEW
376 KNEW = 13 + IZNEW
377
378C REMAINING NUCLEUS GETS A CHARGE WHICH IS ABOUT HALF THE MASS NUMBER
379 ELSEIF ( IANEW .GT. 1 ) THEN
380 IZNEW = FLOAT(IANEW) / 2.15D0 + 0.7D0
381 INACTZ = MAX( IZ - IZNEW, 0 )
382 KNEW = IANEW*100 + IZNEW
383
384C REMAINING NUCLEUS DEEXCITES BY EVAPORATION OF NUCLEONS/ALPHA PARTCLS.
385 IF ( NFRAGM .GE. 2 ) THEN
386 JFIN=0
387 CALL VAPOR(IA,KNEW,JFIN,ITYP,PFRX,PFRY)
388 IF ( JFIN .LE. 0 ) GOTO 190
389 KNEW = 0
390 DO 135 J=1,JFIN
391 EA = GAMMA * PAMA(ITYP(J))
392 IF (DEBUG) WRITE (MDEBUG,*)'SDPM : J,ITYP,EA=',
393 * J,ITYP,SNGL(EA)
394 PTM = EA**2 - PAMA(ITYP(J))**2
395 PT2 = PFRX(J)**2 + PFRY(J)**2
396 IF ( PT2 .GE. PTM ) THEN
397 IF (DEBUG) WRITE(MDEBUG,*)'SDPM : PT REJECT ',J
398 GOTO 135
399 ENDIF
400 IF ( PTM .GT. 0.D0 ) THEN
401 COSTET = SQRT( 1.D0 - PT2/PTM )
402 ELSE
403 COSTET = 1.D0
404 ENDIF
405 IF ( PFRX(J) .NE. 0.D0 .OR. PFRY(J) .NE. 0.D0 ) THEN
406 PHIV = ATAN2( PFRY(J), PFRX(J) )
407 ELSE
408 PHIV = 0.D0
409 ENDIF
410 CALL ADDANG( COSTHE,PHI, COSTET,PHIV,
411 * SECPAR(3),SECPAR(4) )
412 IF ( SECPAR(3) .GE. C(29) ) THEN
413 IF ( J .LT. JFIN ) THEN
414 SECPAR(1) = ITYP(J)
415 CALL TSTACK
416 ELSE
417 KNEW = ITYP(JFIN)
418 IANEW = KNEW/100
419 ENDIF
420 ELSE
421 IF (DEBUG) WRITE(MDEBUG,*)'SDPM : ANGLE REJECT ',J
422 ENDIF
423 135 CONTINUE
424 ENDIF
425 ENDIF
426
427C REMAINING NUCLEUS: MASS 5 CANNOT BE TREATED IN BOX2
428 IF ( KNEW/100 .EQ. 5 ) THEN
429 IF ( MOD(KNEW,100) .GE. 3 ) THEN
430C MASS 5: SPLIT OFF ONE PROTON
431 SECPAR(1) = 14.D0
432 CALL TSTACK
433 KNEW = KNEW - 101
434 ELSE
435C MASS 5: SPLIT OFF ONE NEUTRON
436 SECPAR(1) = 13.D0
437 CALL TSTACK
438 KNEW = KNEW - 100
439 ENDIF
440
441C REMAINING NUCLEUS: MASS 8 CANNOT BE TREATED IN BOX2
442 ELSEIF ( KNEW/100 .EQ. 8 ) THEN
443 IF ( MOD(KNEW,100) .GE. 5 ) THEN
444C MASS 8: SPLIT OFF ONE PROTON
445 SECPAR(1) = 14.D0
446 CALL TSTACK
447 KNEW = KNEW - 101
448 ELSEIF ( MOD(KNEW,100) .LE. 3 ) THEN
449C MASS 8: SPLIT OFF ONE NEUTRON
450 SECPAR(1) = 13.D0
451 CALL TSTACK
452 KNEW = KNEW - 100
453 ELSE
454C MASS 8: SPLIT OFF ONE ALPHA PARTICLE
455 SECPAR(1) = 402.D0
456 CALL TSTACK
457 KNEW = KNEW - 402
458 ENDIF
459 ENDIF
460
461 IF ( KNEW .GT. 0 ) THEN
462 SECPAR(1) = KNEW
463 CALL TSTACK
464 IF ( DEBUG ) WRITE(MDEBUG,555) (SECPAR(I),I=1,9)
465 555 FORMAT (' SDPM : SECPAR=',1P,9E10.3)
466 ENDIF
467 ENDIF
468
469C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
470C HERE THE REACTING NUCLEONS ARE TREATED
471 190 NPRPRO = INACTZ
472 NNEPRO = INACTA - INACTZ
473
474C TREAT INTERACTING NEUTRONS FROM PROJECTILE
475 IF ( NNEPRO .GE. 1 ) THEN
476 CURPAR(1) = 13.D0
477 ITYPE = 13
478C CALCULATE GAMMA, BETA AND ENERGY IN CENTER OF MASS
479 GCM = SQRT( GAMMA * 0.5D0 + 0.5D0 )
480 ECM = PAMA(ITYPE) * GCM * 2.D0
481 BETACM = SQRT( 1.D0 - 1.D0 / GCM**2 )
482 DO 200 LL = 1,NNEPRO
483C USE GHEISHA IF THE CROSS SECTION HAS BEEN CALCULATED FOR GHEISHA
484 IF ( GHEISH .AND. ECM .LE. HILOECM ) THEN
485 ELAB = PAMA(ITYPE) * GAMMA
486 PLAB = ELAB * BETA
487 EKIN = ELAB - PAMA(ITYPE)
488 SIGAIR = CGHSIG(SNGL(PLAB),SNGL(EKIN),ITYPE)
489 CALL CGHEI
490 ELSE
491C DETERMINE TYPE OF INTERACTION FOR NUCLEONS AND ANTINUCLEONS
492 IF ( ECM .GT. CC(4) ) THEN
493C DUAL PARTON MODEL
494 CALL HDPM
495 ELSEIF ( ECM .GT. CC(3) ) THEN
496C USE THE INTERACTION ROUTINES OF PKF GRIEDER
497C 2 HEAVY ISOBARS AND ANNIHILATION
498 CALL BOX63
499 ELSEIF ( ECM .GT. CC(2) ) THEN
500C 1 HEAVY ISOBAR + NUCLEON AND ANNIHILATION
501 CALL BOX62
502 ELSEIF ( ECM .GT. CC(1) ) THEN
503C 1 LIGHT ISOBAR + NUCLEON AND ANNIHILATION
504 CALL BOX61
505 ELSE
506C ELASTIC SCATTERING AND ANNIHILATION
507 CALL BOX60
508 ENDIF
509 ENDIF
510 200 CONTINUE
511 ENDIF
512
513C TREAT INTERACTING PROTONS FROM PROJECTILE IN ROUTINE HDPM
514 IF ( NPRPRO .GE. 1 ) THEN
515 CURPAR(1) = 14.D0
516 ITYPE = 14
517C CALCULATE GAMMA, BETA AND ENERGY IN CENTER OF MASS
518 GCM = SQRT( GAMMA * 0.5D0 + 0.5D0 )
519 ECM = PAMA(ITYPE) * GCM * 2.D0
520 BETACM = SQRT( 1.D0 - 1.D0 / GCM**2 )
521 DO 210 LL = 1,NPRPRO
522C USE GHEISHA IF THE CROSS SECTION HAS BEEN CALCULATED FOR GHEISHA
523 IF ( GHEISH .AND. ECM .LE. HILOECM ) THEN
524 ELAB = PAMA(ITYPE) * GAMMA
525 PLAB = ELAB * BETA
526 EKIN = ELAB - PAMA(ITYPE)
527 SIGAIR = CGHSIG(SNGL(PLAB),SNGL(EKIN),ITYPE)
528 CALL CGHEI
529 ELSE
530C DETERMINE TYPE OF INTERACTION FOR NUCLEONS AND ANTINUCLEONS
531 IF ( ECM .GT. CC(4) ) THEN
532C DUAL PARTON MODEL
533 CALL HDPM
534 ELSEIF ( ECM .GT. CC(3) ) THEN
535C USE THE INTERACTION ROUTINES OF PKF GRIEDER
536C 2 HEAVY ISOBARS AND ANNIHILATION
537 CALL BOX63
538 ELSEIF ( ECM .GT. CC(2) ) THEN
539C 1 HEAVY ISOBAR + NUCLEON AND ANNIHILATION
540 CALL BOX62
541 ELSEIF ( ECM .GT. CC(1) ) THEN
542C 1 LIGHT ISOBAR + NUCLEON AND ANNIHILATION
543 CALL BOX61
544 ELSE
545C ELASTIC SCATTERING AND ANNIHILATION
546 CALL BOX60
547 ENDIF
548 ENDIF
549 210 CONTINUE
550 ENDIF
551
552C ALL PARTICLES, INCLUDING THE LEADING ONE, ARE NOW WRITTEN TO STACK
553
554 ELSE
555 WRITE(MONIOU,*) 'SDPM : NOT FORESEEN PARTICLE TYPE=',ITYPE
556 STOP
557 ENDIF
558
559 RETURN
560 END
Note: See TracBrowser for help on using the repository browser.