COSMOS v7.655  COSMOSv7655
(AirShowerMC)
cvmDecay.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine crhodc (vm, a, np)
 
subroutine comgdc (vm, a, np)
 
subroutine c3pidc (vm, a, np)
 
subroutine cphidc (vm, a, np)
 

Function/Subroutine Documentation

◆ c3pidc()

subroutine c3pidc ( type(ptcl), intent(in)  vm,
type(ptcl), dimension(*), intent(out)  a,
integer, intent(out)  np 
)

Definition at line 95 of file cvmDecay.f.

References cibst1(), cmkptc(), cnbdcy(), and kpion.

Referenced by comgdc(), and cphidc().

95 ! * ****************************
96  implicit none
97 #include "Zptcl.h"
98 #include "Zcode.h"
99  type(ptcl),intent(in):: vm
100  type(ptcl),intent(out):: a(*)
101 ! record /ptcl/ vm, a(*)
102  integer,intent(out):: np
103 
104  integer i
105  integer icon, ntry
106  real(8):: w
107  character*70 msg
108 
109  call cmkptc(kpion, 0, 1, a(1))
110  call cmkptc(kpion, 0, 0, a(2))
111  call cmkptc(kpion, 0, -1, a(3))
112 ! pi+ pi- pi0
113  ntry = 0
114  icon = 1
115  do while( icon .ne. 0 .and. ntry < 200)
116  call cnbdcy(3, vm%mass, a, 0, w, icon)
117  ntry = ntry + 1
118  enddo
119  if(icon .ne. 0) then
120  write(0, *) ' cnbdcy fails for omega 3 pi',
121  * 'cms=',vm%mass, 'icon=',icon
122  np=0
123  else
124  np=3
125  do i = 1, np
126  call cibst1(1, vm, a(i), a(i))
127  enddo
128 !////////////
129 ! write(*,*) ' omega 3body '
130 !//////////////////
131  endif
nodes i
subroutine cibst1(init, p1, p2, po)
Definition: cibst1.f:29
subroutine cnbdcy(n, ecm, p, jw, w, icon)
Definition: cnbdcy.f:48
real(4), save a
Definition: cNRLAtmos.f:20
subroutine cmkptc(code, subcode, charge, p)
Definition: cmkptc.f:15
Definition: Zptcl.h:75
max ptcl codes in the kpion
Definition: Zcode.h:2
Here is the call graph for this function:
Here is the caller graph for this function:

◆ comgdc()

subroutine comgdc ( type(ptcl), intent(in)  vm,
type(ptcl), dimension(*), intent(out)  a,
integer, intent(out)  np 
)

Definition at line 39 of file cvmDecay.f.

References c2bdcy(), c3pidc(), cibst1(), cmkptc(), d0, kmuon, kphoton, kpion, and rndc().

Referenced by cinteomega(), and cvmdcy().

39 ! ******************************
40  implicit none
41 #include "Zptcl.h"
42 #include "Zcode.h"
43 ! record /ptcl/ vm, a(*)
44  type(ptcl),intent(in):: vm
45  type(ptcl),intent(out):: a(*)
46 
47  integer,intent(out):: np
48 
49 
50  real(8):: u
51  integer:: i
52 
53 ! omega--->pi+ pi- pi0 : 88.8 %-->89.2 (2016)
54 ! pi0 gamma: 8.5 -->8.28 (//) 97.48
55 ! pi+ pi- : 2.7 -->1.53 (//) 99.01
56 ! mu+ mu- : 9.0e-5 99.01009
57 ! pi0 mu+ mu-: 1.3e-4 99.01022
58  call rndc(u)
59  if(u .lt. .892d0) then
60  call c3pidc(vm, a, np)
61  elseif(u .lt. .9748d0) then
62  call cmkptc(kpion, 0, 0, a(1))
63  call cmkptc(kphoton, 0, 0, a(2))
64  call c2bdcy(vm, a(1), a(2))
65  np = 2
66  elseif( u < 0.9901d0) then
67 ! pi+ pi-
68  call cmkptc(kpion, 0, 1, a(1))
69  call cmkptc(kpion, 0,-1, a(2))
70  call c2bdcy(vm, a(1), a(2))
71  np = 2
72  elseif( u < 0.9901009d0) then
73 ! mu+ mu-
74  call cmkptc(kmuon, 0, 1, a(1))
75  call cmkptc(kmuon, 0, -1, a(2))
76  call c2bdcy(vm, a(1), a(2))
77  np = 2
78  elseif( u < 0.9901022d0) then
79 ! mu+ mu- pi0
80  call cmkptc(kmuon, 0, 1, a(1))
81  call cmkptc(kmuon, 0, -1, a(2))
82  call cmkptc(kpion, 0, 01, a(3))
83  np = 3
84  do i = 1, np
85  call cibst1(i, vm, a(i), a(i))
86  enddo
87  else
88  ! rest: regards 3pi
89  call c3pidc(vm, a, np)
90  endif
nodes i
const int kphoton
Definition: Zcode.h:6
subroutine cibst1(init, p1, p2, po)
Definition: cibst1.f:29
subroutine rndc(u)
Definition: rnd.f:91
subroutine c3pidc(vm, a, np)
Definition: cvmDecay.f:95
block data cblkEvhnp ! currently usable models data RegMdls ad *special data *Cekaon d0
Definition: cblkEvhnp.h:5
subroutine c2bdcy(p, p1, p2)
Definition: c2bdcy.f:44
real(4), save a
Definition: cNRLAtmos.f:20
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
Here is the call graph for this function:
Here is the caller graph for this function:

◆ cphidc()

subroutine cphidc ( type(ptcl), intent(in)  vm,
type(ptcl), dimension(*), intent(out)  a,
integer, intent(out)  np 
)

Definition at line 134 of file cvmDecay.f.

References c2bdcy(), c3pidc(), cmkptc(), crhodc(), d0, k0l, k0s, keta, kkaon, kmuon, kphoton, kpion, krho, and rndc().

Referenced by cintephi(), and cvmdcy().

134  implicit none
135 #include "Zptcl.h"
136 #include "Zcode.h"
137 
138 ! record /ptcl/ vm, a(*)
139  type(ptcl),intent(in):: vm
140  type(ptcl),intent(out):: a(*)
141 
142  integer,intent(out):: np
143 
144  real(8):: u
145  integer ic, nx
146 ! record /ptcl/pw
147  type(ptcl):: pw
148  nx = 0
149 ! phai-->k+ k- 49.5 %--> 48.9 (2016) %
150 ! k0l k0s 34.4 --> 34.2 (//) 83.1 %
151 ! rho pi 12.9 (rho+ pi-, rho- pi+, rho0 pi0)
152 ! pi+ pi- pi0 1.9 %
153 ! sum of above two 15.32(2016) : so
154 ! we use 12.9-->12.9 1.9->2.42 96.0 %
155 ! 98.42 %
156 ! eta + gamma 1.309 % 99.729
157 ! mu+ mu- 2.87e-4 0.99729
158 ! +0.000287 = 0.997577
159 
160  call rndc(u)
161  if(u .lt. 0.495d0) then
162  call cmkptc(kkaon, 0, 1, a(1))
163  call cmkptc(kkaon, 0, -1, a(2))
164  call c2bdcy(vm, a(1), a(2))
165  np = 2
166  elseif(u .lt. 0.831d0) then
167  call cmkptc(kkaon, k0l, 0, a(1))
168  call cmkptc(kkaon, k0s, 0, a(2))
169  call c2bdcy(vm, a(1), a(2))
170  np = 2
171  elseif(u .lt. 0.960d0) then
172  ! rho + pi
173  call rndc(u)
174  ic=int(3.*u)-1
175  call cmkptc(kpion, 0, ic, a(1))
176  call cmkptc(krho, 0, -ic, a(2))
177  call c2bdcy(vm, a(1), a(2))
178 ! rho is made to decay
179  pw = a(2)
180  call crhodc(pw, a(2), nx)
181  np=nx+1
182  elseif( u < 0.9842d0) then
183  call c3pidc(vm, a, np)
184  elseif( u < 0.99729d0) then
185  ! eta gamma
186  call cmkptc(keta, 0, 0, a(1))
187  call cmkptc(kphoton, 0, 0, a(2))
188  call c2bdcy(vm, a(1), a(2))
189  np =2
190  elseif( u < 0.997577d0 ) then
191  ! mu+ mu-
192  call cmkptc(kmuon, 0, 1, a(1))
193  call cmkptc(kmuon, 0, -1, a(2))
194  call c2bdcy(vm, a(1), a(2))
195  np =2
196  else
197 ! rest is assigned to main br
198  call cmkptc(kkaon, 0, 1, a(1))
199  call cmkptc(kkaon, 0, -1, a(2))
200  call c2bdcy(vm, a(1), a(2))
201  np = 2
202  endif
subroutine crhodc(vm, a, np)
Definition: cvmDecay.f:4
max ptcl codes in the kseethru ! subcode integer k0l
Definition: Zcode.h:2
max ptcl codes in the kseethru ! subcode integer k0s
Definition: Zcode.h:2
const int kphoton
Definition: Zcode.h:6
max ptcl codes in the kkaon
Definition: Zcode.h:2
subroutine rndc(u)
Definition: rnd.f:91
real(8), save pw
Definition: csoftenPiK.f:36
subroutine c3pidc(vm, a, np)
Definition: cvmDecay.f:95
block data cblkEvhnp ! currently usable models data RegMdls ad *special data *Cekaon d0
Definition: cblkEvhnp.h:5
max ptcl codes in the krho
Definition: Zcode.h:2
subroutine c2bdcy(p, p1, p2)
Definition: c2bdcy.f:44
real(4), save a
Definition: cNRLAtmos.f:20
subroutine cmkptc(code, subcode, charge, p)
Definition: cmkptc.f:15
max ptcl codes in the keta
Definition: Zcode.h:2
Definition: Zptcl.h:75
max ptcl codes in the kpion
Definition: Zcode.h:2
max ptcl codes in the kmuon
Definition: Zcode.h:2
Here is the call graph for this function:
Here is the caller graph for this function:

◆ crhodc()

subroutine crhodc ( type(ptcl), intent(in)  vm,
type(ptcl), dimension(*), intent(out)  a,
integer, intent(out)  np 
)

Definition at line 4 of file cvmDecay.f.

References c2bdcy(), cmkptc(), d, kmuon, kpion, and rndc().

Referenced by cinterho(), cphidc(), and cvmdcy().

4 ! ******************************
5  implicit none
6 #include "Zptcl.h"
7 #include "Zcode.h"
8 ! record /ptcl/ vm, a(*)
9  type(ptcl),intent(in):: vm
10  type(ptcl),intent(out):: a(*)
11 
12  integer,intent(out):: np
13 
14  integer:: icharge
15  real(8):: u
16 
17  call rndc(u)
18  if(vm%charge .eq. 0) then
19  if(u < 4.55d-5) then
20 ! mu+ mu- channel
21  call cmkptc(kmuon, 0, 1, a(1))
22  call cmkptc(kmuon, 0, -1, a(2))
23  else
24 ! pi+ pi-
25  call cmkptc(kpion, 0, 1, a(1))
26  call cmkptc(kpion, 0, -1, a(2))
27  endif
28  else
29  ! rho+/-
30  call cmkptc(kpion, 0, 0, a(1))
31  icharge = vm%charge
32  call cmkptc(kpion, 0, icharge, a(2))
33  endif
34  call c2bdcy(vm, a(1), a(2))
35  np = 2
subroutine rndc(u)
Definition: rnd.f:91
dE dx *! Nuc Int sampling table d
Definition: cblkMuInt.h:130
subroutine c2bdcy(p, p1, p2)
Definition: c2bdcy.f:44
real(4), save a
Definition: cNRLAtmos.f:20
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
Here is the call graph for this function:
Here is the caller graph for this function: