COSMOS v7.655  COSMOSv7655
(AirShowerMC)
cetaDecay.f
Go to the documentation of this file.
1 ! testing cetaDecay
2 ! implicit none
3 !#include "Zptcl.h"
4 !#include "Zcode.h"
5 ! record /ptcl/ eta
6 ! record /ptcl/a(10)
7 ! integer n, i, j
8 !c make eta
9 ! call cmkptc(keta, 0, 0, eta)
10 ! eta.fm.p(1) = 1.
11 ! eta.fm.p(2) = -1.
12 ! eta.fm.p(3) = 100.
13 ! eta.fm.p(4) = sqrt(eta.fm.p(3)**2 +
14 ! * eta.fm.p(1)**2 + eta.fm.p(2)**2 + eta.mass**2)
15 ! do i = 1, 10000
16 ! call cetaDecay(eta, a, n)
17 ! if(n .ne. 2) then
18 ! do j =1 , n
19 ! write(*, *) sngl(a(j).fm.p(4)), a(j).code
20 ! enddo
21 ! endif
22 ! enddo
23 ! end
24 ! ******************************************************************
25 ! * *
26 ! * cetaDecay: eta decay
27 ! * *
28 ! ******************************************************************
29 ! eta--> gg 39.41 %
30 ! --> 3 pi0 31.8 % --> 32.68 (2016)
31 ! --> pi+ pi- pi0 23.7 %-->22.92 (2016)
32 ! --> pi+pi- g 4.22% (2016)
33 ! --> mu+ + mu- + gamma 3.1x10^-4 (not %)
34 ! --> mu+ + mu- 5.8 10^-6 (not %)
35 
36  subroutine cetadecay(pj, a, np)
37  implicit none
38 #include "Zptcl.h"
39 #include "Zcode.h"
40 
41  integer np !output. no. of ptcls produced
42 ! record /ptcl/ pj ! input. eta
43  type(ptcl),intent(in):: pj
44 ! record /ptcl/ a(*) ! output. produced ptcls
45  type(ptcl),intent(out):: a(*)
46  integer i, icon
47  real*8 u, w
48 !
49  call rndc(u)
50  if(u .lt. .3941d0) then
51 ! gg
52  call cmkptc(kphoton, 0, 0, a(1))
53  call cmkptc(kphoton, 0, 0, a(2))
54  call c2bdcy(pj, a(1), a(2))
55  np=2
56  elseif(u .lt. .7209d0 ) then
57 ! 3 pi0
58  call cmkptc(kpion, 0, 0, a(1))
59  call cmkptc(kpion, 0, 0, a(2))
60  call cmkptc(kpion, 0, 0, a(3))
61  call cnbdcy(3, pj%mass, a, 0, w, icon)
62  np = 3
63  do i=1, np
64  call cibst1(i, pj, a(i), a(i))
65  enddo
66  elseif( u .lt. 0.9501d0 ) then
67 ! pi+ pi- pi0
68  call cmkptc(kpion, 0, 1, a(1))
69  call cmkptc(kpion, 0, -1, a(2))
70  call cmkptc(kpion, 0, 0, a(3))
71  call cnbdcy(3, pj%mass, a, 0, w, icon)
72  np = 3
73  do i=1, np
74  call cibst1(i, pj, a(i), a(i))
75  enddo
76  elseif( u .lt. 0.9923d0) then
77 ! pi+ pi- + g
78  call cmkptc(kpion, 0, 1, a(1))
79  call cmkptc(kpion, 0, -1, a(2))
80  call cmkptc(kphoton, 0, 0, a(3))
81  call cnbdcy(3, pj%mass, a, 0, w, icon)
82  np = 3
83  do i=1, np
84  call cibst1(i, pj, a(i), a(i))
85  enddo
86  elseif( u .lt. 0.99261d0) then
87 ! mu+ + mu- + g
88  call cmkptc(kmuon, 0, 1, a(1))
89  call cmkptc(kmuon, 0, -1, a(2))
90  call cmkptc(kphoton, 0, 0, a(3))
91  call cnbdcy(3, pj%mass, a, 0, w, icon)
92  np = 3
93  do i=1, np
94  call cibst1(i, pj, a(i), a(i))
95  enddo
96  elseif( u .lt. 0.9926158d0) then
97 ! mu+ + mu-
98  call cmkptc(kmuon, 0, 1, a(1))
99  call cmkptc(kmuon, 0, -1, a(2))
100  np = 2
101  call c2bdcy(pj, a(1), a(2))
102  else
103 ! all rest : regards e+e-g ~ 6.9e-3
104  call cmkptc(kelec, 0, 1, a(1))
105  call cmkptc(kelec, 0, -1, a(2))
106  call cmkptc(kphoton, 0, 0, a(3))
107  call cnbdcy(3, pj%mass, a, 0, w, icon)
108  np = 3
109  do i=1, np
110  call cibst1(i, pj, a(i), a(i))
111  enddo
112  endif
113  end
const int kphoton
Definition: Zcode.h:6
subroutine cibst1(init, p1, p2, po)
Definition: cibst1.f:29
max ptcl codes in the kelec
Definition: Zcode.h:2
subroutine cnbdcy(n, ecm, p, jw, w, icon)
Definition: cnbdcy.f:48
subroutine rndc(u)
Definition: rnd.f:91
block data cblkEvhnp ! currently usable models data RegMdls ad *special data *Cekaon d0
Definition: cblkEvhnp.h:5
subroutine cetadecay(pj, a, np)
Definition: cetaDecay.f:37
subroutine c2bdcy(p, p1, p2)
Definition: c2bdcy.f:44
subroutine cmkptc(code, subcode, charge, p)
Definition: cmkptc.f:15
Definition: Zptcl.h:75
max ptcl codes in the kpion
Definition: Zcode.h:2
max ptcl codes in the kmuon
Definition: Zcode.h:2