COSMOS v7.655  COSMOSv7655
(AirShowerMC)
cdDecay.f
Go to the documentation of this file.
1 ! ******************************************************************
2 ! * *
3 ! * prompt muon decay mode is important.
4 ! * others are roughly introduced to keep the energy conservation
5 ! *
6 !
7 ! ******************************************************************
8 ! ******** note: c2bdcy: boost to lab is done inside
9 ! cnbdcy: boost must be done outside.
10  subroutine cddecay(pj, a, np)
11  implicit none
12 #include "Zptcl.h"
13 #include "Zcode.h"
14 
15  type(ptcl):: pj ! input. demeson
16  integer np ! output. # of ptcls stored in a
17  type(ptcl):: a(*) ! output. to store produced ptcls
18 
19  real*8 u, sumbr
20 
21 
22  call rndc(u)
23  if(pj%charge .ne. 0) then
24  if(u .lt. 0.093d0 ) then
25 ! K0b + mu + neu 9.3%
26  call cddecay1( pj, a, np)
27  elseif(u .lt. 0.132d0) then
28 ! K + pi + mu + neumu 3.9 % 13.2
29  call cddecay2(pj, a, np)
30  else
31  sumbr=0.4179 ! the total % below)/100.
32 ! K 2pi 9.22 9.22
33 ! K0b + e + neu 8.6 % 17.82
34 ! K0s + 2pi 6.8 24.62
35 ! K+3pi 6.0 30.62
36 ! K + pi+ e+ neue 4.1 % 34.72
37 ! K0s + 3pi 3.02 37.74
38 ! K0L + pi 1.46 39.20
39 ! K0s + pi 1.45 40.65
40 ! 3pi + pi0 1.14 41.79
41  call rndc(u)
42 ! u=u*sumbr ! 2013. enu prob too big
43  if(u .lt. 0.0922) then
44 ! D+ K- pi+ p+ 9.22 9.22
45  call cddecay3(pj, a, np)
46  elseif(u .lt. 0.1782) then
47 ! K0b e+ neu 8.6 % 17.82
48  call cddecay4(pj, a, np)
49  elseif(u .lt. 0.2462) then
50 ! K0s + 2pi 6.8 24.62
51  call cddecay5( pj, a, np)
52  elseif(u .lt. 0.3062) then
53 ! K+3pi 6.0 30.62
54  call cddecay6(pj, a, np)
55  elseif(u .lt. 0.3472) then
56 ! K + pi+ e+ neue 4.1 % 34.72
57  call cddecay7( pj, a, np)
58  elseif( u .lt. 0.3774) then
59 ! K0s + 3pi 3.02 37.74
60  call cddecay8(pj, a, np)
61  elseif( u .lt. 0.3920) then
62 ! K0L + pi 1.46 39.20
63  call cddecay9(pj, a, np)
64  elseif( u .lt. 0.4065) then
65 ! K0s + pi 1.45 40.65
66  call cddecay10(pj, a, np)
67  elseif( u < 0.4179) then
68 ! 3pi + pi0 1.14 41.79
69  call cddecay11(pj, a, np)
70  else ! 2013.Mar 2. neglect other
71  np =0
72  endif
73  endif
74  else
75 ! D0
76 ! K^- + mu^+ + neu 3.31 0.0331
77 ! K^- + e^+ + neue 3.55 0.06863
78 ! K*+ e^+ + neue 2.16
79 ! treat as ==> K+ pi + e^+ + neue and combine with
80 ! K^-+ pi0+e^+ +neue 1.6
81 ! K^0bar + pi^-+e^+ +neue 2.7 0.1332
82 !c K*+ mu^+ + neu 1.91
83 ! treat as K+pi+ mu^+ + neu 0.1523
84 ! -----------------------
85 ! K- pi+ pi- 13.9 0.2913
86 ! K-+2pi+ pi- 8.08 0.3721
87 ! K-2pi+ pi- pi0 4.2 0.4141
88 ! K- rho+=>K- pi+ pi0 10.8 0.5221
89 
90 ! K0s + pi+ pi- pi0 5.2 0.5741
91 
92 ! ............
93 
94  call rndc(u)
95  if(u .lt. 0.0331) then
96 ! D0 K- mu+ neu 3.31
97  call cddecay20(pj, a, np)
98  elseif( u < 0.06863) then
99 ! K + e + neue 3.55
100  call cddecay25( pj, a, np)
101  elseif( u < 0.1332) then
102  call cddecay28(pj, a, np)
103  elseif( u< 0.2523) then
104 ! treat as K+pi+ mu^+ + neu 0.1523
105  call cddecay29(pj, a, np)
106  else
107 ! call rndc(u)
108 ! sumbr =0.3695
109 ! u=u*sumbr ! avoid too large other branch
110  if(u .lt. 0.2913) then
111 ! K- pi+ pi- 13.9 13.9
112  call cddecay21(pj, a, np)
113  elseif(u .lt. 0.3721 ) then
114 ! K-+2pi+ pi- 8.08
115  call cddecay22( pj, a, np)
116  elseif(u .lt. 0.4141) then
117 ! K-2pi+ pi- pi0 4.2 23.5
118  call cddecay23( pj, a, np)
119  elseif(u .lt. 0.5221) then
120 ! K- pi+ pi0 10.8 27.39
121  call cddecay24( pj, a, np)
122  elseif(u .lt. 0.5741 ) then
123 ! K0s + pi+ pi- pi0 5.2
124  call cddecay27(pj, a, np)
125  else ! neglect all others
126  np = 0
127  endif
128  endif
129  endif
130  end
131 ! ********************************
132  subroutine cddecay1( pj, a, np)
133  implicit none
134 #include "Zptcl.h"
135 #include "Zcode.h"
136 ! D+-> K0b + mu + neu 9.3%
137  type(ptcl):: pj ! input. demeson
138  integer np ! output. # of ptcls stored in a
139  type(ptcl):: a(*) ! output. to store produced ptcls
140 
141  integer muchg, nusubc, kchg, ksubc, icon, i, echg
142  integer pichg
143 
144  real*8 u, w
145 
146  call rndc(u)
147  if(u .lt. .50) then
148  ksubc = k0s
149  else
150  ksubc = k0l
151  endif
152  muchg = pj%charge
153  nusubc = muchg
154 ! neue
155  call cmkptc(kneue, nusubc, 0, a(1))
156  call cmkptc(kmuon, 0, muchg, a(2))
157  call cmkptc(kkaon, ksubc, 0, a(3))
158 ! 3 body pure phase space
159  call cnbdcy(3, pj%mass, a, 0, w, icon)
160  np = 3
161  do i=1, np
162  call cibst1(i, pj, a(i), a(i))
163  enddo
164  return
165 ! *******************************************
166 ! D+ K- pi+ mu+ neumu 3.9 % 13.2
167  entry cddecay2(pj, a, np)
168 
169  kchg =-pj%charge
170  nusubc = regptcl
171  call cmkptc(kkaon, 0, kchg, a(1))
172  call cmkptc(kpion, 0, -kchg, a(2))
173  call cmkptc(kmuon, 0, -kchg, a(3))
174  call cmkptc(kneumu, nusubc, 0, a(4))
175 
176  call cnbdcy(4, pj%mass, a, 0, w, icon)
177  np = 4
178  do i=1, np
179  call cibst1(i, pj, a(i), a(i))
180  enddo
181  return
182 ! ***********************************
183 ! D+ K- pi+ p+ 9.22 9.22
184  entry cddecay3(pj, a, np)
185 
186 
187  kchg = -pj%charge
188  call cmkptc(kkaon, 0, kchg, a(1))
189  call cmkptc(kpion, 0, -kchg, a(2))
190  call cmkptc(kpion, 0, -kchg, a(3))
191 
192  call cnbdcy(3, pj%mass, a, 0, w, icon)
193  np = 3
194  do i=1, np
195  call cibst1(i, pj, a(i), a(i))
196  enddo
197  return
198 !*****************************************
199 ! K0b e+ neu 8.6 % 17.82
200  entry cddecay4(pj, a, np)
201 
202 
203  ksubc = k0s
204  echg = 1
205  call rndc(u)
206  if(u .lt.0.5 ) then
207  ksubc = k0l
208  echg = -1
209  endif
210  nusubc = -echg
211  call cmkptc(kkaon, ksubc, 0, a(1))
212  call cmkptc(kelec, 0, echg, a(2))
213  call cmkptc(kneue, nusubc, 0, a(3))
214 
215  call cnbdcy(3, pj%mass, a, 0, w, icon)
216  np = 3
217  do i=1, np
218  call cibst1(i, pj, a(i), a(i))
219  enddo
220  return
221 ! *******************************
222 ! K0s pi+ pi0 6.8 24.62
223  entry cddecay5( pj, a, np)
224 
225 
226  ksubc = k0s
227  pichg = pj%charge
228  call cmkptc(kkaon, ksubc, 0, a(1))
229  call cmkptc(kpion, 0, pichg, a(2))
230  call cmkptc(kpion, 0, 0, a(3))
231 
232  call cnbdcy(3, pj%mass, a, 0, w, icon)
233  np = 3
234  do i=1, np
235  call cibst1(i, pj, a(i), a(i))
236  enddo
237  return
238 ! *************************************
239 ! K-pi+ pi+ pi0 6.0 30.62
240  entry cddecay6( pj, a, np)
241 
242  kchg = -pj%charge
243  call cmkptc(kkaon, 0, kchg, a(1))
244  call cmkptc(kpion, 0, -kchg, a(2))
245  call cmkptc(kpion, 0, -kchg, a(3))
246  call cmkptc(kpion, 0, 0, a(4))
247 
248  call cnbdcy(4, pj%mass, a, 0, w, icon)
249  np = 4
250  do i=1, np
251  call cibst1(i, pj, a(i), a(i))
252  enddo
253  return
254 ! **********************************
255 ! K- pi+ e+ neue 4.1 % 34.72
256  entry cddecay7( pj, a, np)
257 
258  kchg =-pj%charge
259  nusubc = -kchg
260  call cmkptc(kkaon, 0, kchg, a(1))
261  call cmkptc(kpion, 0, -kchg, a(2))
262  call cmkptc(kelec, 0, -kchg, a(3))
263  call cmkptc(kneue, nusubc, 0, a(4))
264 
265  call cnbdcy(4, pj%mass, a, 0, w, icon)
266  np = 4
267  do i=1, np
268  call cibst1(i, pj, a(i), a(i))
269  enddo
270  return
271 !*****************************
272 ! K0s + pi+ pi+ pi- 3.02 37.74
273  entry cddecay8(pj, a, np)
274 
275  ksubc = k0s
276  pichg = pj%charge
277  call cmkptc(kkaon, ksubc, 0, a(1))
278  call cmkptc(kpion, 0, pichg, a(2))
279  call cmkptc(kpion, 0, pichg, a(3))
280  call cmkptc(kpion, 0, -pichg, a(4))
281 
282 
283  call cnbdcy(4, pj%mass, a, 0, w, icon)
284  np = 4
285  do i=1, np
286  call cibst1(i, pj, a(i), a(i))
287  enddo
288  return
289 ! ******************************
290 ! K0L + pi 1.46 39.20
291  entry cddecay9(pj, a, np)
292 
293  ksubc = k0l
294  pichg = pj%charge
295  call cmkptc(kkaon, ksubc, 0, a(1))
296  call cmkptc(kpion, 0, pichg, a(2))
297  call c2bdcy(pj, a(1), a(2))
298  np =2
299  return
300 ! ****************************
301 ! K0s + pi 1.46 39.20
302  entry cddecay10(pj, a, np)
303 
304  ksubc = k0s
305  pichg = pj%charge
306  call cmkptc(kkaon, ksubc, 0, a(1))
307  call cmkptc(kpion, 0, pichg, a(2))
308  call c2bdcy(pj, a(1), a(2))
309  np =2
310  return
311 ! ***************************
312 ! 3pi + pi0 1.14 41.79
313  entry cddecay11(pj, a, np)
314 
315  pichg = pj%charge
316  call cmkptc(kpion, 0, pichg, a(1))
317  call cmkptc(kpion, 0, pichg, a(2))
318  call cmkptc(kpion, 0, -pichg, a(3))
319  call cmkptc(kpion, 0, 0, a(4))
320 
321  call cnbdcy(4, pj%mass, a, 0, w, icon)
322  np = 4
323  do i=1, np
324  call cibst1(i, pj, a(i), a(i))
325  enddo
326  return
327 ! **************************************
328 ! D0 K- mu+ neu 3.31
329  entry cddecay20(pj, a, np)
330 
331 
332  call cmkptc(kkaon, 0, -1, a(1))
333  call cmkptc(kmuon, 1, 1, a(2))
334  call cmkptc(kneumu, -1, 0, a(3))
335 
336  call cnbdcy(3, pj%mass, a, 0, w, icon)
337  np = 3
338  do i=1, np
339  call cibst1(i, pj, a(i), a(i))
340  enddo
341  return
342 !********************************************
343 ! D0 K- pi+ pi- 13.9 13.9
344  entry cddecay21(pj, a, np)
345 
346  call cmkptc(kkaon, 0, -1, a(1))
347  call cmkptc(kpion, 0, 1, a(2))
348  call cmkptc(kpion, 0, -1, a(3))
349 
350  call cnbdcy(3, pj%mass, a, 0, w, icon)
351  np = 3
352  do i=1, np
353  call cibst1(i, pj, a(i), a(i))
354  enddo
355  return
356 !****************************************************
357 ! K-+2pi+ pi- 5.4 19.3
358  entry cddecay22( pj, a, np)
359 
360 
361  call cmkptc(kkaon, 0, -1, a(1))
362  call cmkptc(kpion, 0, 1, a(2))
363  call cmkptc(kpion, 0, 1, a(3))
364  call cmkptc(kpion, 0, -1, a(4))
365 
366  call cnbdcy(4, pj%mass, a, 0, w, icon)
367  np = 4
368  do i=1, np
369  call cibst1(i, pj, a(i), a(i))
370  enddo
371  return
372 ! *********************************
373 ! K-2pi+ pi- pi0 4.2 23.5
374  entry cddecay23( pj, a, np)
375 
376 
377  call cmkptc(kkaon, 0, -1, a(1))
378  call cmkptc(kpion, 0, 1, a(2))
379  call cmkptc(kpion, 0, 1, a(3))
380  call cmkptc(kpion, 0, -1, a(4))
381  call cmkptc(kpion, 0, 0, a(5))
382 
383  call cnbdcy(5, pj%mass, a, 0, w, icon)
384  np = 5
385  do i=1, np
386  call cibst1(i, pj, a(i), a(i))
387  enddo
388  return
389 ! ***************************
390 ! K- pi+ pi0 3.89 27.39
391  entry cddecay24( pj, a, np)
392 
393  call cmkptc(kkaon, 0, -1, a(1))
394  call cmkptc(kpion, 0, 1, a(2))
395  call cmkptc(kpion, 0, 0, a(3))
396  call cnbdcy(3, pj%mass, a, 0, w, icon)
397  np =3
398  do i=1, np
399  call cibst1(i, pj, a(i), a(i))
400  enddo
401 
402  return
403 ! *************************
404 ! K + e + neue 3.55
405  entry cddecay25( pj, a, np)
406  if( pj%subcode == antip ) then
407  kchg = -1
408  echg = 1
409  else
410  kchg = 1
411  echg = -1
412  endif
413  call cmkptc(kkaon, kchg, kchg, a(1))
414  call cmkptc(kelec, echg, echg, a(2))
415  call cmkptc(kneue, kchg, 0, a(3))
416 
417  call cnbdcy(3, pj%mass, a, 0, w, icon)
418  np = 3
419  do i=1, np
420  call cibst1(i, pj, a(i), a(i))
421  enddo
422  return
423 !**************************************
424 ! K0s + pi+ pi- 2.99 33.96
425  entry cddecay26(pj, a, np)
426 
427 
428  call cmkptc(kkaon, 0, -1, a(1))
429  call cmkptc(kpion, 0, 1, a(2))
430  call cmkptc(kpion, 0, -1, a(3))
431 
432  call cnbdcy(3, pj%mass, a, 0, w, icon)
433  np = 3
434  do i=1, np
435  call cibst1(i, pj, a(i), a(i))
436  enddo
437  return
438 ! *******************************
439 ! K0s + pi+ pi- pi0 5.2
440  entry cddecay27(pj, a, np)
441  ksubc = k0s
442  call cmkptc(kkaon, ksubc, 0, a(1))
443  call cmkptc(kpion, 0, 1, a(2))
444  call cmkptc(kpion, 0, -1, a(3))
445  call cmkptc(kpion, 0, 0, a(4))
446 
447  call cnbdcy(4, pj%mass, a, 0, w, icon)
448  np = 4
449  do i=1, np
450  call cibst1(i, pj, a(i), a(i))
451  enddo
452 ! *******************************
453 ! K + pi + e^+ +neue
454  entry cddecay28(pj, a, np)
455  call rndc(u)
456  if( pj%subcode == antip ) then
457  echg = 1
458  else
459  echg = -1
460  endif
461  if( u < 0.5 ) then
462  kchg = 0
463  pichg = -echg
464  else
465  kchg = -echg
466  pichg = 0
467  endif
468  call cmkptc(kkaon, 0, kchg, a(1))
469  call cmkptc(kpion, 0, pichg, a(2))
470  call cmkptc(kelec, 0, echg, a(3))
471  call cmkptc(kneue, -echg, 0, a(4))
472 
473  call cnbdcy(4, pj%mass, a, 0, w, icon)
474  np = 4
475  do i=1, np
476  call cibst1(i, pj, a(i), a(i))
477  enddo
478  return
479 
480 ! *******************************
481 ! K+pi+ mu^+ + neu 0.1523
482  entry cddecay29(pj, a, np)
483  call rndc(u)
484  if( pj%subcode == antip ) then
485  muchg = 1
486  else
487  muchg = -1
488  endif
489  if( u < 0.5 ) then
490  kchg = 0
491  pichg = -muchg
492  else
493  kchg = -muchg
494  pichg = 0
495  endif
496  call cmkptc(kkaon, 0, kchg, a(1))
497  call cmkptc(kpion, 0, pichg, a(2))
498  call cmkptc(kmuon, 0, muchg, a(3))
499  call cmkptc(kneumu, -muchg, 0, a(4))
500 
501  call cnbdcy(4, pj%mass, a, 0, w, icon)
502  np = 4
503  do i=1, np
504  call cibst1(i, pj, a(i), a(i))
505  enddo
506  end
cddecay
subroutine cddecay(pj, a, np)
Definition: cdDecay.f:11
cddecay1
subroutine cddecay1(pj, a, np)
Definition: cdDecay.f:133
k0l
max ptcl codes in the kseethru ! subcode integer k0l
Definition: Zcode.h:2
k0s
max ptcl codes in the kseethru ! subcode integer k0s
Definition: Zcode.h:2
kkaon
max ptcl codes in the kkaon
Definition: Zcode.h:2
cibst1
subroutine cibst1(init, p1, p2, po)
Definition: cibst1.f:29
kelec
max ptcl codes in the kelec
Definition: Zcode.h:2
kneue
max ptcl codes in the kneue
Definition: Zcode.h:2
cnbdcy
subroutine cnbdcy(n, ecm, p, jw, w, icon)
Definition: cnbdcy.f:48
rndc
subroutine rndc(u)
Definition: rnd.f:91
regptcl
max ptcl codes in the kseethru ! subcode integer regptcl
Definition: Zcode.h:2
d0
block data cblkEvhnp ! currently usable models data RegMdls ad *special data *Cekaon d0
Definition: cblkEvhnp.h:5
kneumu
max ptcl codes in the kneumu
Definition: Zcode.h:2
c2bdcy
subroutine c2bdcy(p, p1, p2)
Definition: c2bdcy.f:44
cmkptc
subroutine cmkptc(code, subcode, charge, p)
Definition: cmkptc.f:15
ptcl
Definition: Zptcl.h:75
antip
max ptcl codes in the kseethru ! subcode integer antip
Definition: Zcode.h:2
kpion
max ptcl codes in the kpion
Definition: Zcode.h:2
kmuon
max ptcl codes in the kmuon
Definition: Zcode.h:2