1 /* copyright(c) 1991-2012 free software foundation, inc.
2 this file is part of the gnu c library.
4 the gnu c library is free software; you can redistribute it and/or
5 modify it under the terms of the gnu lesser general
Public
6 license as published by the free software foundation; either
7 version 2.1 of the license, or(at your option) any later version.
9 the gnu c library is distributed in the hope that it will be useful,
10 but without any warranty; without even the implied warranty of
11 merchantability or fitness for a particular purpose. see the gnu
12 lesser general
Public license for more details.
14 you should have received a copy of the gnu lesser general
Public
15 license along with the gnu c library;
if not, see
16 <http://www.gnu.org/licenses/>. */
19 /* this header is separate from features.h so that the compiler can
20 include it implicitly at the start of every compilation. it must
21 not itself include <features.h> or any other header that includes
22 <features.h> because the
implicit include comes before any feature
23 test macros that may be defined in a source file before it first
24 explicitly includes a system header. gcc knows the name of this
25 header in order to preinclude it. */
27 /* we
do support the iec 559 math functionality,
real and complex. */
29 /* wchar_t uses iso/iec 10646 (2nd ed., published 2011-03-15) /
32 /* we
do not support c11 <threads.h>. */
34 */////////////////////////////////////////////////////////////////////////////////////
38 *// due to short
common block names it owerwrites variables in other parts //
41 *// one should add suffix c_taul_ to names of all commons as soon as possible
43 */////////////////////////////////////////////////////////////////////////////////////
45 */////////////////////////////////////////////////////////////////////////////////////
47 *// standard tauola interface/initialization routines of functionality exactly //
48 *// as in tauola cpc but input is partially from xpar(*) matrix //
49 *// itauxpar is for indirect adressing //
51 */////////////////////////////////////////////////////////////////////////////////////
54 SUBROUTINE inietc(ITAUXPAR,xpar)
62 COMMON / taurad / xk0dec,itdkrc
63 DOUBLE PRECISION xk0dec
64 COMMON / jaki / jak1,jak2,jakp,jakm,ktom
65 * note: i dont see keya1=2,3 realy implemented in the code sj. ??????
69 * keya1=1 constant width of a1 and rho
70 * keya1=2 free choice of rho propagator(defined in
function fpik)
71 * and free choice of a1 mass and width.
function g(Q**2)
72 * (see formula 3.48 in comp. phys. comm. 64 (1991) 275)
73 * hard coded both in monte carlo and in testing distribution.
74 * keya1=3
function g(Q**2) hardcoded in the Monte Carlo
75 * (it is timy to calculate
77 idff = xpar(itauxpar+3)
79 xk0dec = xpar(itauxpar+5)
80 c radiative correction switch in tau --> e(mu) decays
81 itdkrc = xpar(itauxpar+4)
82 c switches of tau+ tau- decay modes
83 jak1 = xpar(itauxpar+1)
84 jak2 = xpar(itauxpar+2)
85 c output file number for tauola
87 c keya1 is used for formfactors actually not in use
88 keya1 = xpar(itauxpar+6)
91 WRITE(iout,bxtxt)
' Parameters passed from KK to Tauola: '
92 WRITE(iout,bxl1i) jak1,
'dec. type 1-st tau ',
'Jak1 ',
't01'
93 WRITE(iout,bxl1i) jak2,
'dec. type 2-nd tau ',
'Jak2 ',
't02'
94 WRITE(iout,bxl1i) keya1,
'current type a1 dec.',
'KeyA1 ',
't03'
95 WRITE(iout,bxl1i) idff,
'PDG id 1-st tau ',
'idff ',
't04'
96 WRITE(iout,bxl1i) itdkrc,
'R.c. switch lept dec',
'itdkRC',
't05'
97 WRITE(iout,bxl1g) xk0dec,
'IR-cut for lept r.c.',
'xk0dec',
't06'
102 SUBROUTINE initdk(ITAUXPAR,xpar)
103 * ----------------------------------------------------------------------
104 * initialisation of tau decay parameters and routines
107 * ----------------------------------------------------------------------
115 COMMON / decpar / gfermi,gv,ga,ccabib,scabib,gamel
116 REAL*4 gfermi,gv,ga,ccabib,scabib,gamel
117 COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
118 * ,ampiz,ampi,amro,gamro,ama1,gama1
119 * ,amk,amkz,amkst,gamkst
121 REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
122 * ,ampiz,ampi,amro,gamro,ama1,gama1
123 * ,amk,amkz,amkst,gamkst
124 COMMON / taubra / gamprt(30),jlist(30),nchan
125 COMMON / taukle / bra1,brk0,brk0b,brks
126 REAL*4 bra1,brk0,brk0b,brks
127 parameter(nmode=15,nm1=0,nm2=1,nm3=8,nm4=2,nm5=1,nm6=3)
128 COMMON / taudcd /idffin(9,nmode),mulpik(nmode)
130 CHARACTER names(nmode)*31
131 CHARACTER oldnames(7)*31
134 $ bxinit =
'(1x,1h*,g17.8, 16x, a31,a4,a4, 1x,1h*)'
139 * list of branching ratios
140 cam normalised to e nu nutau channel
141 cam enu munu pinu rhonu a1nu knu k*nu pi
142 cam
DATA jlist / 1, 2, 3, 4, 5, 6, 7,
143 *am
DATA gamprt /1.000,0.9730,0.6054,1.2432,0.8432,0.0432,o.o811,0.616
147 * conventions of particles names
148 * k-,p-,k+, k0,p-,kb, k-,p0,k0
149 * 3, 1,-3 , 4, 1,-4 , 3, 2, 4 ,
150 * p0,p0,k-, k-,p-,p+, p-,kb,p0
151 * 2, 2, 3 , 3, 1,-1 , 1,-4, 2 ,
156 dimension nopik(6,nmode),npik(nmode)
157 *am outgoing multiplicity and flavors of multi-pion /multi-k modes
166 DATA nopik / -1,-1, 1, 2, 0, 0, 2, 2, 2,-1, 0, 0,
167 1 -1,-1, 1, 2, 2, 0, -1,-1,-1, 1, 1, 0,
168 2 -1,-1,-1, 1, 1, 2, -1,-1, 1, 2, 2, 2,
169 3 -3,-1, 3, 0, 0, 0, -4,-1, 4, 0, 0, 0,
170 4 -3, 2,-4, 0, 0, 0, 2, 2,-3, 0, 0, 0,
171 5 -3,-1, 1, 0, 0, 0, -1, 4, 2, 0, 0, 0,
172 6 9,-1, 2, 0, 0, 0, -1, 2, 8, 0, 0, 0,
173 c ajwmod fix sign bug, 2/22/99
174 7 -3,-4, 0, 0, 0, 0 /
175 * list of branching ratios
180 IF(i.EQ. 1) gamprt(i) =0.1800
181 IF(i.EQ. 2) gamprt(i) =0.1751
182 IF(i.EQ. 3) gamprt(i) =0.1110
183 IF(i.EQ. 4) gamprt(i) =0.2515
184 IF(i.EQ. 5) gamprt(i) =0.1790
185 IF(i.EQ. 6) gamprt(i) =0.0071
186 IF(i.EQ. 7) gamprt(i) =0.0134
187 IF(i.EQ. 8) gamprt(i) =0.0450
188 IF(i.EQ. 9) gamprt(i) =0.0100
189 IF(i.EQ.10) gamprt(i) =0.0009
190 IF(i.EQ.11) gamprt(i) =0.0004
191 IF(i.EQ.12) gamprt(i) =0.0003
192 IF(i.EQ.13) gamprt(i) =0.0005
193 IF(i.EQ.14) gamprt(i) =0.0015
194 IF(i.EQ.15) gamprt(i) =0.0015
195 IF(i.EQ.16) gamprt(i) =0.0015
196 IF(i.EQ.17) gamprt(i) =0.0005
197 IF(i.EQ.18) gamprt(i) =0.0050
198 IF(i.EQ.19) gamprt(i) =0.0055
199 IF(i.EQ.20) gamprt(i) =0.0017
200 IF(i.EQ.21) gamprt(i) =0.0013
201 IF(i.EQ.22) gamprt(i) =0.0010
202 IF(i.EQ. 1) oldnames(i)=
' TAU- --> E- '
203 IF(i.EQ. 2) oldnames(i)=
' TAU- --> MU- '
204 IF(i.EQ. 3) oldnames(i)=
' TAU- --> PI- '
205 IF(i.EQ. 4) oldnames(i)=
' TAU- --> PI-, PI0 '
206 IF(i.EQ. 5) oldnames(i)=
' TAU- --> A1- (two subch) '
207 IF(i.EQ. 6) oldnames(i)=
' TAU- --> K- '
208 IF(i.EQ. 7) oldnames(i)=
' TAU- --> K*- (two subch) '
209 IF(i.EQ. 8) names(i-7)=
' TAU- --> 2PI-, PI0, PI+ '
210 IF(i.EQ. 9) names(i-7)=
' TAU- --> 3PI0, PI- '
211 IF(i.EQ.10) names(i-7)=
' TAU- --> 2PI-, PI+, 2PI0 '
212 IF(i.EQ.11) names(i-7)=
' TAU- --> 3PI-, 2PI+, '
213 IF(i.EQ.12) names(i-7)=
' TAU- --> 3PI-, 2PI+, PI0 '
214 IF(i.EQ.13) names(i-7)=
' TAU- --> 2PI-, PI+, 3PI0 '
215 IF(i.EQ.14) names(i-7)=
' TAU- --> K-, PI-, K+ '
216 IF(i.EQ.15) names(i-7)=
' TAU- --> K0, PI-, K0B '
217 IF(i.EQ.16) names(i-7)=
' TAU- --> K-, K0, PI0 '
218 IF(i.EQ.17) names(i-7)=
' TAU- --> PI0 PI0 K- '
219 IF(i.EQ.18) names(i-7)=
' TAU- --> K- PI- PI+ '
220 IF(i.EQ.19) names(i-7)=
' TAU- --> PI- K0B PI0 '
221 IF(i.EQ.20) names(i-7)=
' TAU- --> ETA PI- PI0 '
222 IF(i.EQ.21) names(i-7)=
' TAU- --> PI- PI0 GAM '
223 IF(i.EQ.22) names(i-7)=
' TAU- --> K- K0 '
232 idffin(j,i)=nopik(j,i)
237 * --- coefficients to fix ratio of:
238 * --- a1 3charged/ a1 1charged 2 neutrals matrix elements(masless lim.)
239 * --- probability of k0 to be ks
240 * --- probability of k0b to be ks
241 * --- ratio of coefficients for k*--> k0 pi-
242 * --- all coefficents should be in the range(0.0,1.0)
243 * --- they meaning is probability of the first choice only
IF one
244 * --- neglects mass-phase space effects
259 IF (xpar(itauxpar+100+1).GT.-1d0)
THEN
260 c initialization form kk
261 ccabib = xpar(itauxpar+7)
262 gv = xpar(itauxpar+8)
263 ga = xpar(itauxpar+9)
265 bra1 = xpar(itauxpar+10)
266 brks = xpar(itauxpar+11)
267 brk0 = xpar(itauxpar+12)
268 brk0b = xpar(itauxpar+13)
270 gamprt(k)=xpar(itauxpar+100+k)
273 * zw 13.04.89 here was an error
274 scabib = sqrt(1.-ccabib**2)
276 gamel = gfermi**2*amtau**5/(192*pi**3)
278 * CALL dexay(-1,pol1)
280 * printouts for kk version
289 WRITE(iout,bxtxt)
' TAUOLA Initialization SUBROUTINE INITDK: '
290 WRITE(iout,bxtxt)
' Adopted to read from KK '
291 WRITE(iout,bxtxt)
' '
292 WRITE(iout,bxtxt)
' Choice Probability -- Decay Channel'
294 WRITE(iout,bxinit) gamprt(k)/sum, oldnames(k),
'****',
'***'
297 WRITE(iout,bxinit) gamprt(k)/sum, names(k-7),
'****',
'***'
299 WRITE(iout,bxtxt)
' In addition:'
300 WRITE(iout,bxinit) gv,
'Vector W-tau-nu coupl. ',
'****',
'***'
301 WRITE(iout,bxinit) ga,
'Axial W-tau-nu coupl. ',
'****',
'***'
302 WRITE(iout,bxinit) gfermi,
'Fermi Coupling ',
'****',
'***'
303 WRITE(iout,bxinit) ccabib,
'cabibo angle ',
'****',
'***'
304 WRITE(iout,bxinit) bra1,
'a1 br ratio (massless) ',
'****',
'***'
305 WRITE(iout,bxinit) brks,
'K* br ratio (massless) ',
'****',
'***'
311 SUBROUTINE iniphy(XK00)
312 * ----------------------------------------------------------------------
313 * initialisation of parameters
314 * used in qed and/or gsw routines
315 * ----------------------------------------------------------------------
316 COMMON / qedprm /alfinv,alfpi,xk0
317 REAL*8 alfinv,alfpi,xk0
320 pi8 = 4.d0*datan(1.d0)
322 alfpi = 1d0/(alfinv*pi8)
326 SUBROUTINE inimas(ITAUXPAR,xpar)
327 * ----------------------------------------------------------------------
328 * initialisation of masses
331 * ----------------------------------------------------------------------
338 COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
339 * ,ampiz,ampi,amro,gamro,ama1,gama1
340 * ,amk,amkz,amkst,gamkst
342 REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
343 * ,ampiz,ampi,amro,gamro,ama1,gama1
344 * ,amk,amkz,amkst,gamkst
347 $ bxinit =
'(1x,1h*,g17.8, 16x, a31,a4,a4, 1x,1h*)'
350 * in-coming / out-going fermion masses
358 * masses used in tau decays
372 c in-coming / out-going fermion masses
377 c masses used in tau decays cleo settings
391 WRITE(iout,bxtxt)
' TAUOLA Initialization SUBROUTINE INIMAS: '
392 WRITE(iout,bxtxt)
' Adopted to read from KK '
393 WRITE(iout,bxinit) amtau,
'AMTAU tau-mass ',
'****',
'***'
394 WRITE(iout,bxinit) amel ,
'AMEL electron-mass ',
'****',
'***'
395 WRITE(iout,bxinit) ammu ,
'AMMU muon-mass ',
'****',
'***'
399 SUBROUTINE choice(MNUM,RR,ICHAN,PROB1,PROB2,PROB3,
400 $ amrx,gamrx,amra,gamra,amrb,gamrb)
401 COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
402 * ,ampiz,ampi,amro,gamro,ama1,gama1
403 * ,amk,amkz,amkst,gamkst
405 REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
406 * ,ampiz,ampi,amro,gamro,ama1,gama1
407 * ,amk,amkz,amkst,gamkst
413 c xxxxa correspond to s2 channel
423 ELSEIF(mnum.EQ.1)
THEN
432 ELSEIF(mnum.EQ.2)
THEN
441 ELSEIF(mnum.EQ.3)
THEN
450 ELSEIF(mnum.EQ.4)
THEN
459 ELSEIF(mnum.EQ.5)
THEN
468 ELSEIF(mnum.EQ.6)
THEN
477 ELSEIF(mnum.EQ.7)
THEN
486 ELSEIF(mnum.EQ.8)
THEN
495 ELSEIF(mnum.EQ.101)
THEN
504 ELSEIF(mnum.EQ.102)
THEN
524 IF (rr.LE.prob1)
THEN
526 ELSEIF(rr.LE.(prob1+prob2))
THEN
541 prob3=1.0-prob1-prob2
543 FUNCTION dcdmas(IDENT)
544 COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
545 * ,ampiz,ampi,amro,gamro,ama1,gama1
546 * ,amk,amkz,amkst,gamkst
548 REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
549 * ,ampiz,ampi,amro,gamro,ama1,gama1
550 * ,amk,amkz,amkst,gamkst
551 IF (ident.EQ. 1)
THEN
553 ELSEIF (ident.EQ.-1)
THEN
555 ELSEIF (ident.EQ. 2)
THEN
557 ELSEIF (ident.EQ.-2)
THEN
559 ELSEIF (ident.EQ. 3)
THEN
561 ELSEIF (ident.EQ.-3)
THEN
563 ELSEIF (ident.EQ. 4)
THEN
565 ELSEIF (ident.EQ.-4)
THEN
567 ELSEIF (ident.EQ. 8)
THEN
569 ELSEIF (ident.EQ.-8)
THEN
571 ELSEIF (ident.EQ. 9)
THEN
573 ELSEIF (ident.EQ.-9)
THEN
576 print *,
'STOP IN APKMAS, WRONG IDENT=',ident
581 FUNCTION lunpik(ID,ISGN)
582 COMMON / taukle / bra1,brk0,brk0b,brks
583 REAL*4 bra1,brk0,brk0b,brks
586 IF (ident.EQ. 1)
THEN
588 ELSEIF (ident.EQ.-1)
THEN
590 ELSEIF (ident.EQ. 2)
THEN
592 ELSEIF (ident.EQ.-2)
THEN
594 ELSEIF (ident.EQ. 3)
THEN
596 ELSEIF (ident.EQ.-3)
THEN
598 ELSEIF (ident.EQ. 4)
THEN
600 * k0 --> k0_long(is 130) / k0_short(is 310) = 1/1
602 IF (xio(1).GT.brk0)
THEN
607 ELSEIF (ident.EQ.-4)
THEN
609 * k0b--> k0_long(is 130) / k0_short(is 310) = 1/1
611 IF (xio(1).GT.brk0b)
THEN
616 ELSEIF (ident.EQ. 8)
THEN
618 ELSEIF (ident.EQ.-8)
THEN
620 ELSEIF (ident.EQ. 9)
THEN
622 ELSEIF (ident.EQ.-9)
THEN
625 print *,
'STOP IN IPKDEF, WRONG IDENT=',ident
633 SUBROUTINE taurdf(KTO)
634 c this routine can be called before any tau+ or tau- event is generated
635 c it can be used to generate tau+ and tau- samples of different
637 COMMON / taukle / bra1,brk0,brk0b,brks
638 REAL*4 bra1,brk0,brk0b,brks
639 COMMON / taubra / gamprt(30),jlist(30),nchan
642 c ajwmod: set the brs for(a1+ -> rho+ pi0) and(k*+ -> k0 pi+)
649 c ajwmod: set the brs for(a1+ -> rho+ pi0) and(k*+ -> k0 pi+)