source: trunk/MagicSoft/Simulation/Corsika/Mmcs/pi0dec.f@ 6525

Last change on this file since 6525 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: 4.5 KB
Line 
1 SUBROUTINE PI0DEC
2
3C-----------------------------------------------------------------------
4C PI 0 DEC(AY)
5C
6C DECAY OF PI0 INTO 2 GAMMAS OR INTO E(+) + E(-) + GAMMA
7C THIS SUBROUTINE IS CALLED FROM BOX3
8C
9C DESIGN : D. HECK IK3 FZK KARLSRUHE
10C-----------------------------------------------------------------------
11
12 IMPLICIT NONE
13*KEEP,CONST.
14 COMMON /CONST/ PI,PI2,OB3,TB3,ENEPER
15 DOUBLE PRECISION PI,PI2,OB3,TB3,ENEPER
16*KEEP,DECAY.
17 COMMON /DECAY/ GAM345,COS345,PHI345
18 DOUBLE PRECISION GAM345(3),COS345(3),PHI345(3)
19*KEEP,GENER.
20 COMMON /GENER/ GEN,ALEVEL
21 DOUBLE PRECISION GEN,ALEVEL
22*KEEP,PAM.
23 COMMON /PAM/ PAMA,SIGNUM
24 DOUBLE PRECISION PAMA(6000),SIGNUM(6000)
25*KEEP,PARPAR.
26 COMMON /PARPAR/ CURPAR,SECPAR,PRMPAR,OUTPAR,C,
27 * E00,E00PN,PTOT0,PTOT0N,THICKH,ITYPE,LEVL
28 DOUBLE PRECISION CURPAR(14),SECPAR(14),PRMPAR(14),OUTPAR(14),
29 * C(50),E00,E00PN,PTOT0,PTOT0N,THICKH
30 INTEGER ITYPE,LEVL
31*KEEP,PARPAE.
32 DOUBLE PRECISION GAMMA,COSTHE,PHI,H,T,X,Y,CHI,BETA,GCM,ECM
33 EQUIVALENCE (CURPAR(2),GAMMA), (CURPAR(3),COSTHE),
34 * (CURPAR(4), PHI ), (CURPAR(5), H ),
35 * (CURPAR(6), T ), (CURPAR(7), X ),
36 * (CURPAR(8), Y ), (CURPAR(9), CHI ),
37 * (CURPAR(10),BETA), (CURPAR(11),GCM ),
38 * (CURPAR(12),ECM )
39*KEEP,RANDPA.
40 COMMON /RANDPA/ FAC,U1,U2,RD,NSEQ,ISEED,KNOR
41 DOUBLE PRECISION FAC,U1,U2
42 REAL RD(3000)
43 INTEGER ISEED(103,10),NSEQ
44 LOGICAL KNOR
45*KEEP,RUNPAR.
46 COMMON /RUNPAR/ FIXHEI,THICK0,HILOECM,HILOELB,
47 * STEPFC,NRRUN,NSHOW,PATAPE,MONIIN,
48 * MONIOU,MDEBUG,NUCNUC,
49 * CETAPE,
50 * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
51 * N1STTR,MDBASE,
52 * DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
53 * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
54 * ,GHEISH,GHESIG
55 COMMON /RUNPAC/ DSN,HOST,USER
56 DOUBLE PRECISION FIXHEI,THICK0,HILOECM,HILOELB
57 REAL STEPFC
58 INTEGER NRRUN,NSHOW,PATAPE,MONIIN,MONIOU,MDEBUG,NUCNUC,
59 * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
60 * N1STTR,MDBASE
61 INTEGER CETAPE
62 CHARACTER*79 DSN
63 CHARACTER*20 HOST,USER
64
65 LOGICAL DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
66 * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
67 * ,GHEISH,GHESIG
68*KEND.
69
70 DOUBLE PRECISION AUX1,AUX2,COSTH1,COSTH2,EPITO2,FI1
71 INTEGER I
72C-----------------------------------------------------------------------
73
74 IF ( DEBUG ) WRITE(MDEBUG,444) (CURPAR(I),I=1,9)
75 444 FORMAT (' PI0DEC: CURPAR=',1P,9E10.3)
76
77C COPY VERTEX COORDINATES INTO SECPAR
78 DO 1 I = 5,8
79 SECPAR(I) = CURPAR(I)
80 1 CONTINUE
81 SECPAR( 9) = GEN
82 SECPAR(10) = ALEVEL
83
84C LOOK FOR DECAY MODE
85 CALL RMMAR (RD,3,1)
86
87C DECAY PI(0) ----> GAMMA + GAMMA
88 IF ( RD(3) .LT. 0.98802 ) THEN
89C HALF OF TOTAL ENERGY OF THE PION = EPITO2
90 EPITO2 = 0.5D0 * GAMMA * PAMA(7)
91 AUX1 = 1.D0 + BETA * RD(1)
92 AUX2 = 1.D0 - BETA * RD(1)
93 COSTH1 = (BETA + RD(1)) / AUX1
94 COSTH2 = (BETA - RD(1)) / AUX2
95
96C FIRST GAMMA (WITH HIGHER ENERGY)
97 FI1 = PI2 * RD(2)
98 CALL ADDANG( COSTHE,PHI, COSTH1,FI1, SECPAR(3),SECPAR(4) )
99 IF ( SECPAR(3) .GT. C(29) ) THEN
100 SECPAR(1) = 1.D0
101C ENERGY OF GAMMA
102 SECPAR(2) = AUX1 * EPITO2
103 CALL TSTACK
104 ENDIF
105
106C SECOND GAMMA (WITH LOWER ENERGY)
107 CALL ADDANG( COSTHE,PHI, COSTH2,FI1+PI, SECPAR(3),SECPAR(4) )
108 IF ( SECPAR(3) .GT. C(29) ) THEN
109 SECPAR(1) = 1.D0
110C ENERGY OF GAMMA
111 SECPAR(2) = AUX2 * EPITO2
112 CALL TSTACK
113 ENDIF
114
115C DECAY PI(0) ----> E(+) + E(-) + GAMMA (DALITZ DECAY)
116C (UNIFORM PHASE SPACE DISTRIBUTION IS ASSUMED FOR THIS DECAY)
117 ELSE
118 CALL DECAY6( PAMA(7), PAMA(3), PAMA(2), 0.D0,
119 * 0.D0,0.D0,0.D0, 1.D0, 2)
120 DO 11 I = 1,3
121 CALL ADDANG( COSTHE,PHI, COS345(I),PHI345(I),
122 * SECPAR(3),SECPAR(4) )
123 IF ( SECPAR(3) .GT. C(29) ) THEN
124 SECPAR(1) = FLOAT(4 - I)
125 SECPAR(2) = GAM345(I)
126 CALL TSTACK
127 ENDIF
128 11 CONTINUE
129
130 ENDIF
131
132 RETURN
133 END
Note: See TracBrowser for help on using the repository browser.