C++InterfacetoTauola
tauola-F/jetset-F/tauola_photos_ini.F
1 C this file is created by hand from taumain.F
2 C actions: Remove routines: TAUDEM DECTES TAUFIL FILHEP
3 C add: INIETC will not necesarily work fine ...
4 C replace TRALO4
5 C rename INIPHY to INIPHX
6 
7  SUBROUTINE inietc(jakk1,jakk2,itd,ifpho)
8  COMMON / idfc / idff
9  COMMON / taurad / xk0dec,itdkrc
10  DOUBLE PRECISION xk0dec
11  COMMON / jaki / jak1,jak2,jakp,jakm,ktom
12  COMMON /phoact/ ifphot
13  SAVE
14 C KTO=1 will denote tau+, thus :: IDFF=-15
15  idff=-15
16 C XK0 for tau decays.
17  xk0dec=0.01
18 C radiative correction switch in tau --> e (mu) decays !
19  itdkrc=itd
20 C switches of tau+ tau- decay modes !!
21  jak1=jakk1
22  jak2=jakk2
23 C photos activation switch
24  ifphot=ifpho
25  end
26 
27  SUBROUTINE tralo4(KTOS,PHOI,PHOF,AM)
28 !! Corrected 11.10.96 (ZW) tralor for KORALW.
29 !! better treatment is to cascade from tau rest-frame through W
30 !! restframe down to LAB.
31  COMMON / momdec / q1,q2,p1,p2,p3,p4
32  COMMON /tralid/ idtra
33  double precision q1(4),q2(4),p1(4),p2(4),p3(4),p4(4),p1qq(4),p2qq(4)
34  double precision pin(4),pout(4),pbst(4),pbs1(4),qq(4),pi
35  double precision thet,phi,exe
36  REAL*4 phoi(4),phof(4)
37  SAVE
38  DATA pi /3.141592653589793238462643d0/
39  am=sqrt(abs
40  $ (phoi(4)**2-phoi(3)**2-phoi(2)**2-phoi(1)**2))
41  idtra=ktos
42  DO k=1,4
43  pin(k)=phoi(k)
44  phof(k)=phoi(k)
45  ENDDO
46 ! write(*,*) idtra
47  IF (idtra.EQ.1) THEN
48  DO k=1,4
49  pbst(k)=p1(k)
50  qq(k)=q1(k)
51  ENDDO
52  ELSEIF(idtra.EQ.2) THEN
53  DO k=1,4
54  pbst(k)=p2(k)
55  qq(k)=q1(k)
56  ENDDO
57  ELSEIF(idtra.EQ.3) THEN
58  DO k=1,4
59  pbst(k)=p3(k)
60  qq(k)=q2(k)
61  ENDDO
62  ELSE
63  DO k=1,4
64  pbst(k)=p4(k)
65  qq(k)=q2(k)
66  ENDDO
67  ENDIF
68 
69 
70 
71  CALL bostdq(1,qq,pbst,pbst)
72  CALL bostdq(1,qq,p1,p1qq)
73  CALL bostdq(1,qq,p2,p2qq)
74  pbs1(4)=pbst(4)
75  pbs1(3)=sqrt(pbst(3)**2+pbst(2)**2+pbst(1)**2)
76  pbs1(2)=0d0
77  pbs1(1)=0d0
78  exe=(pbs1(4)+pbs1(3))/dsqrt(pbs1(4)**2-pbs1(3)**2)
79 C for KTOS=1 boost is antiparallel to 4-momentum of P2.
80 C restframes of tau+ tau- and 'first' frame of 'higgs' are all connected
81 C by boosts along z axis
82  IF(ktos.EQ.1) exe=(pbs1(4)-pbs1(3))/dsqrt(pbs1(4)**2-pbs1(3)**2)
83  CALL bostd3(exe,pin,pout)
84 
85 C once in Z/gamma/Higgs rest frame we control further kinematics by P2QQ for KTOS=1,2
86  thet=acos(p2qq(3)/sqrt(p2qq(3)**2+p2qq(2)**2+p2qq(1)**2))
87  phi=0d0
88  phi=acos(p2qq(1)/sqrt(p2qq(2)**2+p2qq(1)**2))
89  IF(p2qq(2).LT.0d0) phi=2*pi-phi
90 
91  CALL rotpox(thet,phi,pout)
92  CALL bostdq(-1,qq,pout,pout)
93  DO k=1,4
94  phof(k)=pout(k)
95  ENDDO
96  END
97 
98 
99  SUBROUTINE choice(MNUM,RR,ICHAN,PROB1,PROB2,PROB3,
100  $ amrx,gamrx,amra,gamra,amrb,gamrb)
101  COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
102  * ,ampiz,ampi,amro,gamro,ama1,gama1
103  * ,amk,amkz,amkst,gamkst
104 C
105  REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
106  * ,ampiz,ampi,amro,gamro,ama1,gama1
107  * ,amk,amkz,amkst,gamkst
108 C
109  amrop=1.1
110  gamrop=0.36
111  amom=.782
112  gamom=0.0084
113 C XXXXA CORRESPOND TO S2 CHANNEL !
114  IF(mnum.EQ.0) THEN
115  prob1=0.5
116  prob2=0.5
117  amrx =ama1
118  gamrx=gama1
119  amra =amro
120  gamra=gamro
121  amrb =amro
122  gamrb=gamro
123  ELSEIF(mnum.EQ.1) THEN
124  prob1=0.5
125  prob2=0.5
126  amrx =1.57
127  gamrx=0.9
128  amrb =amkst
129  gamrb=gamkst
130  amra =amro
131  gamra=gamro
132  ELSEIF(mnum.EQ.2) THEN
133  prob1=0.5
134  prob2=0.5
135  amrx =1.57
136  gamrx=0.9
137  amrb =amkst
138  gamrb=gamkst
139  amra =amro
140  gamra=gamro
141  ELSEIF(mnum.EQ.3) THEN
142  prob1=0.5
143  prob2=0.5
144  amrx =1.27
145  gamrx=0.3
146  amra =amkst
147  gamra=gamkst
148  amrb =amkst
149  gamrb=gamkst
150  ELSEIF(mnum.EQ.4) THEN
151  prob1=0.5
152  prob2=0.5
153  amrx =1.27
154  gamrx=0.3
155  amra =amkst
156  gamra=gamkst
157  amrb =amkst
158  gamrb=gamkst
159  ELSEIF(mnum.EQ.5) THEN
160  prob1=0.5
161  prob2=0.5
162  amrx =1.27
163  gamrx=0.3
164  amra =amkst
165  gamra=gamkst
166  amrb =amro
167  gamrb=gamro
168  ELSEIF(mnum.EQ.6) THEN
169  prob1=0.4
170  prob2=0.4
171  amrx =1.27
172  gamrx=0.3
173  amra =amro
174  gamra=gamro
175  amrb =amkst
176  gamrb=gamkst
177  ELSEIF(mnum.EQ.7) THEN
178  prob1=0.0
179  prob2=1.0
180  amrx =1.27
181  gamrx=0.9
182  amra =amro
183  gamra=gamro
184  amrb =amro
185  gamrb=gamro
186  ELSEIF(mnum.EQ.8) THEN
187  prob1=0.0
188  prob2=1.0
189  amrx =amrop
190  gamrx=gamrop
191  amrb =amom
192  gamrb=gamom
193  amra =amro
194  gamra=gamro
195  ELSEIF(mnum.EQ.101) THEN
196  prob1=.35
197  prob2=.35
198  amrx =1.2
199  gamrx=.46
200  amrb =amom
201  gamrb=gamom
202  amra =amom
203  gamra=gamom
204  ELSEIF(mnum.EQ.102) THEN
205  prob1=0.0
206  prob2=0.0
207  amrx =1.4
208  gamrx=.6
209  amrb =amom
210  gamrb=gamom
211  amra =amom
212  gamra=gamom
213  ELSE
214  prob1=0.0
215  prob2=0.0
216  amrx =ama1
217  gamrx=gama1
218  amra =amro
219  gamra=gamro
220  amrb =amro
221  gamrb=gamro
222  ENDIF
223 C
224  IF (rr.LE.prob1) THEN
225  ichan=1
226  ELSEIF(rr.LE.(prob1+prob2)) THEN
227  ichan=2
228  ax =amra
229  gx =gamra
230  amra =amrb
231  gamra=gamrb
232  amrb =ax
233  gamrb=gx
234  px =prob1
235  prob1=prob2
236  prob2=px
237  ELSE
238  ichan=3
239  ENDIF
240 C
241  prob3=1.0-prob1-prob2
242  END
243  SUBROUTINE initdk
244 * ----------------------------------------------------------------------
245 * INITIALISATION OF TAU DECAY PARAMETERS and routines
246 *
247 * called by : KORALZ
248 * ----------------------------------------------------------------------
249 
250  COMMON / decpar / gfermi,gv,ga,ccabib,scabib,gamel
251  REAL*4 gfermi,gv,ga,ccabib,scabib,gamel
252  COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
253  * ,ampiz,ampi,amro,gamro,ama1,gama1
254  * ,amk,amkz,amkst,gamkst
255 *
256  REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
257  * ,ampiz,ampi,amro,gamro,ama1,gama1
258  * ,amk,amkz,amkst,gamkst
259  COMMON / taubra / gamprt(30),jlist(30),nchan
260  COMMON / taukle / bra1,brk0,brk0b,brks
261  REAL*4 bra1,brk0,brk0b,brks
262 #if defined (ALEPH)
263  parameter(nmode=15,nm1=0,nm2=1,nm3=8,nm4=2,nm5=1,nm6=3)
264  COMMON / taudcd /idffin(9,nmode),mulpik(nmode)
265  & ,names
266  CHARACTER names(nmode)*31
267 #else
268  parameter(nmode=15,nm1=0,nm2=1,nm3=8,nm4=2,nm5=1,nm6=3)
269  COMMON / taudcd /idffin(9,nmode),mulpik(nmode)
270  & ,names
271  CHARACTER names(nmode)*31
272 #endif
273  CHARACTER oldnames(7)*31
274  CHARACTER*80 bxinit
275  parameter(
276  $ bxinit ='(1x,1h*,g17.8, 16x, a31,a4,a4, 1x,1h*)'
277  $ )
278  REAL*4 pi,pol1(4)
279 *
280 *
281 * LIST OF BRANCHING RATIOS
282 CAM normalised to e nu nutau channel
283 CAM enu munu pinu rhonu A1nu Knu K*nu pi
284 CAM DATA JLIST / 1, 2, 3, 4, 5, 6, 7,
285 #if defined (ALEPH)
286 CAM /0.1779,0.1731,0.1106,0.2530,0.1811,0.0072,0.0139
287 CAM DATA GAMPRT / 1.000,0.9732,0.6217,1.4221,1.0180,0.0405,0.0781
288 CAM DATA GAMPRT /1.000,0.9676,0.6154,1.3503,1.0225,0.0368,O.O758
289 CAM
290 C
291 C conventions of particles names
292 c
293 cam mode (JAK) 8 9
294 CAM channel pi- pi- pi0 pi+ 3pi0 pi-
295 cam particle code -1,-1, 2, 1, 0, 0, 2, 2, 2,-1, 0, 0,
296 CAM BR relative to electron .2414, .0601,
297 c
298 * 10 11
299 * 1 3pi+- 2pi0 5pi+-
300 * 1 -1,-1, 1, 2, 2, 0, -1,-1,-1, 1, 1, 0,
301 * 1 .0281, .0045,
302 
303 * 12 13
304 * 2 5pi+- pi0 3pi+- 3pi0
305 * 2 -1,-1,-1, 1, 1, 2, -1,-1, 1, 2, 2, 2,
306 * 2 .0010, .0062,
307 
308 * 14 15
309 * 3 K- pi- K+ K0 pi- KB
310 * 3 -3,-1, 3, 0, 0, 0, 4,-1,-4, 0, 0, 0,
311 * 3 .0096, .0169,
312 
313 * 16 17
314 * 4 K- pi0 K0 2pi0 K-
315 * 4 -3, 2, 4, 0, 0, 0, 2, 2,-3, 0, 0, 0,
316 * 4 .0056, .0045,
317 
318 * 18 19
319 * 5 K- pi- pi+ pi- KB pi0
320 * 5 -3,-1, 1, 0, 0, 0, -1,-4, 2, 0, 0, 0,
321 * 5 .0219, .0180,
322 
323 * 20 21
324 * 6 eta pi- pi0 pi- pi0 gamma
325 * 6 9,-1, 2, 0, 0, 0, -1, 2, 8, 0, 0, 0,
326 * 6 .0096, .0088,
327 
328 * 22 /
329 * 7 K- K0 /
330 * 7 -3, 4 /
331 * 7 .0146 /
332 #else
333 *AM DATA GAMPRT /1.000,0.9730,0.6054,1.2432,0.8432,0.0432,O.O811,0.616
334 *AM
335 *AM multipion decays
336 *
337 * conventions of particles names
338 * K-,P-,K+, K0,P-,KB, K-,P0,K0
339 * 3, 1,-3 , 4, 1,-4 , 3, 2, 4 ,
340 * P0,P0,K-, K-,P-,P+, P-,KB,P0
341 * 2, 2, 3 , 3, 1,-1 , 1,-4, 2 ,
342 * ET,P-,P0 P-,P0,GM
343 * 9, 1, 2 , 1, 2, 8
344 *
345 #endif
346 C
347  dimension nopik(6,nmode),npik(nmode)
348 *AM outgoing multiplicity and flavors of multi-pion /multi-K modes
349  DATA npik / 4, 4,
350  1 5, 5,
351  2 6, 6,
352  3 3, 3,
353  4 3, 3,
354  5 3, 3,
355  6 3, 3,
356  7 2 /
357 #if defined (ALEPH)
358  DATA nopik / -1,-1, 2, 1, 0, 0, 2, 2, 2,-1, 0, 0,
359  1 -1,-1, 1, 2, 2, 0, -1,-1,-1, 1, 1, 0,
360  2 -1,-1,-1, 1, 1, 2, -1,-1, 1, 2, 2, 2,
361  3 -3,-1, 3, 0, 0, 0, 4,-1,-4, 0, 0, 0,
362  4 -3, 2, 4, 0, 0, 0, 2, 2,-3, 0, 0, 0,
363  5 -3,-1, 1, 0, 0, 0, -1,-4, 2, 0, 0, 0,
364  6 9,-1, 2, 0, 0, 0, -1, 2, 8, 0, 0, 0,
365 #else
366  DATA nopik / -1,-1, 1, 2, 0, 0, 2, 2, 2,-1, 0, 0,
367  1 -1,-1, 1, 2, 2, 0, -1,-1,-1, 1, 1, 0,
368  2 -1,-1,-1, 1, 1, 2, -1,-1, 1, 2, 2, 2,
369  3 -3,-1, 3, 0, 0, 0, -4,-1, 4, 0, 0, 0,
370  4 -3, 2,-4, 0, 0, 0, 2, 2,-3, 0, 0, 0,
371  5 -3,-1, 1, 0, 0, 0, -1, 4, 2, 0, 0, 0,
372  6 9,-1, 2, 0, 0, 0, -1, 2, 8, 0, 0, 0,
373 #endif
374 #if defined (CLEO)
375 C AJWMOD fix sign bug, 2/22/99
376  7 -3,-4, 0, 0, 0, 0 /
377 #else
378  7 -3, 4, 0, 0, 0, 0 /
379 #endif
380 * LIST OF BRANCHING RATIOS
381  nchan = nmode + 7
382  DO 1 i = 1,30
383  IF (i.LE.nchan) THEN
384  jlist(i) = i
385 #if defined (CePeCe)
386  IF(i.EQ. 1) gamprt(i) = 1.0000
387  IF(i.EQ. 2) gamprt(i) = 1.0000
388  IF(i.EQ. 3) gamprt(i) = 1.0000
389  IF(i.EQ. 4) gamprt(i) = 1.0000
390  IF(i.EQ. 5) gamprt(i) = 1.0000
391  IF(i.EQ. 6) gamprt(i) = 1.0000
392  IF(i.EQ. 7) gamprt(i) = 1.0000
393  IF(i.EQ. 8) gamprt(i) = 1.0000
394  IF(i.EQ. 9) gamprt(i) = 1.0000
395  IF(i.EQ.10) gamprt(i) = 1.0000
396  IF(i.EQ.11) gamprt(i) = 1.0000
397  IF(i.EQ.12) gamprt(i) = 1.0000
398  IF(i.EQ.13) gamprt(i) = 1.0000
399  IF(i.EQ.14) gamprt(i) = 1.0000
400  IF(i.EQ.15) gamprt(i) = 1.0000
401  IF(i.EQ.16) gamprt(i) = 1.0000
402  IF(i.EQ.17) gamprt(i) = 1.0000
403  IF(i.EQ.18) gamprt(i) = 1.0000
404  IF(i.EQ.19) gamprt(i) = 1.0000
405  IF(i.EQ.20) gamprt(i) = 1.0000
406  IF(i.EQ.21) gamprt(i) = 1.0000
407  IF(i.EQ.22) gamprt(i) = 1.0000
408 #elif defined (CLEO)
409  IF(i.EQ. 1) gamprt(i) =0.1800
410  IF(i.EQ. 2) gamprt(i) =0.1751
411  IF(i.EQ. 3) gamprt(i) =0.1110
412  IF(i.EQ. 4) gamprt(i) =0.2515
413  IF(i.EQ. 5) gamprt(i) =0.1790
414  IF(i.EQ. 6) gamprt(i) =0.0071
415  IF(i.EQ. 7) gamprt(i) =0.0134
416  IF(i.EQ. 8) gamprt(i) =0.0450
417  IF(i.EQ. 9) gamprt(i) =0.0100
418  IF(i.EQ.10) gamprt(i) =0.0009
419  IF(i.EQ.11) gamprt(i) =0.0004
420  IF(i.EQ.12) gamprt(i) =0.0003
421  IF(i.EQ.13) gamprt(i) =0.0005
422  IF(i.EQ.14) gamprt(i) =0.0015
423  IF(i.EQ.15) gamprt(i) =0.0015
424  IF(i.EQ.16) gamprt(i) =0.0015
425  IF(i.EQ.17) gamprt(i) =0.0005
426  IF(i.EQ.18) gamprt(i) =0.0050
427  IF(i.EQ.19) gamprt(i) =0.0055
428  IF(i.EQ.20) gamprt(i) =0.0017
429  IF(i.EQ.21) gamprt(i) =0.0013
430  IF(i.EQ.22) gamprt(i) =0.0010
431 #elif defined (ALEPH)
432  IF(i.EQ. 1) gamprt(i) = 1.0000
433  IF(i.EQ. 2) gamprt(i) = .9732
434  IF(i.EQ. 3) gamprt(i) = .6217
435  IF(i.EQ. 4) gamprt(i) = 1.4221
436  IF(i.EQ. 5) gamprt(i) = 1.0180
437  IF(i.EQ. 6) gamprt(i) = .0405
438  IF(i.EQ. 7) gamprt(i) = .0781
439  IF(i.EQ. 8) gamprt(i) = .2414
440  IF(i.EQ. 9) gamprt(i) = .0601
441  IF(i.EQ.10) gamprt(i) = .0281
442  IF(i.EQ.11) gamprt(i) = .0045
443  IF(i.EQ.12) gamprt(i) = .0010
444  IF(i.EQ.13) gamprt(i) = .0062
445  IF(i.EQ.14) gamprt(i) = .0096
446  IF(i.EQ.15) gamprt(i) = .0169
447  IF(i.EQ.16) gamprt(i) = .0056
448  IF(i.EQ.17) gamprt(i) = .0045
449  IF(i.EQ.18) gamprt(i) = .0219
450  IF(i.EQ.19) gamprt(i) = .0180
451  IF(i.EQ.20) gamprt(i) = .0096
452  IF(i.EQ.21) gamprt(i) = .0088
453  IF(i.EQ.22) gamprt(i) = .0146
454 #else
455 #endif
456  IF(i.EQ. 1) oldnames(i)=' TAU- --> E- '
457  IF(i.EQ. 2) oldnames(i)=' TAU- --> MU- '
458  IF(i.EQ. 3) oldnames(i)=' TAU- --> PI- '
459  IF(i.EQ. 4) oldnames(i)=' TAU- --> PI-, PI0 '
460  IF(i.EQ. 5) oldnames(i)=' TAU- --> A1- (two subch) '
461  IF(i.EQ. 6) oldnames(i)=' TAU- --> K- '
462  IF(i.EQ. 7) oldnames(i)=' TAU- --> K*- (two subch) '
463  IF(i.EQ. 8) names(i-7)=' TAU- --> 2PI-, PI0, PI+ '
464  IF(i.EQ. 9) names(i-7)=' TAU- --> 3PI0, PI- '
465  IF(i.EQ.10) names(i-7)=' TAU- --> 2PI-, PI+, 2PI0 '
466  IF(i.EQ.11) names(i-7)=' TAU- --> 3PI-, 2PI+, '
467  IF(i.EQ.12) names(i-7)=' TAU- --> 3PI-, 2PI+, PI0 '
468  IF(i.EQ.13) names(i-7)=' TAU- --> 2PI-, PI+, 3PI0 '
469  IF(i.EQ.14) names(i-7)=' TAU- --> K-, PI-, K+ '
470  IF(i.EQ.15) names(i-7)=' TAU- --> K0, PI-, K0B '
471 #if defined (ALEPH)
472  IF(i.EQ.16) names(i-7)=' TAU- --> K- PI0 K0 '
473 #else
474  IF(i.EQ.16) names(i-7)=' TAU- --> K-, K0, PI0 '
475 #endif
476  IF(i.EQ.17) names(i-7)=' TAU- --> PI0 PI0 K- '
477  IF(i.EQ.18) names(i-7)=' TAU- --> K- PI- PI+ '
478  IF(i.EQ.19) names(i-7)=' TAU- --> PI- K0B PI0 '
479  IF(i.EQ.20) names(i-7)=' TAU- --> ETA PI- PI0 '
480  IF(i.EQ.21) names(i-7)=' TAU- --> PI- PI0 GAM '
481  IF(i.EQ.22) names(i-7)=' TAU- --> K- K0 '
482  ELSE
483  jlist(i) = 0
484  gamprt(i) = 0.
485  ENDIF
486  1 CONTINUE
487  DO i=1,nmode
488  mulpik(i)=npik(i)
489  DO j=1,mulpik(i)
490  idffin(j,i)=nopik(j,i)
491  ENDDO
492  ENDDO
493 *
494 *
495 * --- COEFFICIENTS TO FIX RATIO OF:
496 * --- A1 3CHARGED/ A1 1CHARGED 2 NEUTRALS MATRIX ELEMENTS (MASLESS LIM.)
497 * --- PROBABILITY OF K0 TO BE KS
498 * --- PROBABILITY OF K0B TO BE KS
499 * --- RATIO OF COEFFICIENTS FOR K*--> K0 PI-
500 * --- ALL COEFFICENTS SHOULD BE IN THE RANGE (0.0,1.0)
501 * --- THEY MEANING IS PROBABILITY OF THE FIRST CHOICE ONLY IF ONE
502 * --- NEGLECTS MASS-PHASE SPACE EFFECTS
503  bra1=0.5
504  brk0=0.5
505  brk0b=0.5
506  brks=0.6667
507 *
508 
509  gfermi = 1.16637e-5
510  ccabib = 0.975
511  gv = 1.0
512  ga =-1.0
513 
514 
515 
516 * ZW 13.04.89 HERE WAS AN ERROR
517  scabib = sqrt(1.-ccabib**2)
518  pi =4.*atan(1.)
519  gamel = gfermi**2*amtau**5/(192*pi**3)
520 *
521  CALL dexay(-1,pol1)
522 *
523  RETURN
524  END
525  FUNCTION dcdmas(IDENT)
526  COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
527  * ,ampiz,ampi,amro,gamro,ama1,gama1
528  * ,amk,amkz,amkst,gamkst
529 *
530  REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
531  * ,ampiz,ampi,amro,gamro,ama1,gama1
532  * ,amk,amkz,amkst,gamkst
533  IF (ident.EQ. 1) THEN
534  apkmas=ampi
535  ELSEIF (ident.EQ.-1) THEN
536  apkmas=ampi
537  ELSEIF (ident.EQ. 2) THEN
538  apkmas=ampiz
539  ELSEIF (ident.EQ.-2) THEN
540  apkmas=ampiz
541  ELSEIF (ident.EQ. 3) THEN
542  apkmas=amk
543  ELSEIF (ident.EQ.-3) THEN
544  apkmas=amk
545  ELSEIF (ident.EQ. 4) THEN
546  apkmas=amkz
547  ELSEIF (ident.EQ.-4) THEN
548  apkmas=amkz
549  ELSEIF (ident.EQ. 8) THEN
550  apkmas=0.0001
551  ELSEIF (ident.EQ.-8) THEN
552  apkmas=0.0001
553  ELSEIF (ident.EQ. 9) THEN
554  apkmas=0.5488
555  ELSEIF (ident.EQ.-9) THEN
556  apkmas=0.5488
557  ELSE
558  print *, 'STOP IN APKMAS, WRONG IDENT=',ident
559  stop
560  ENDIF
561  dcdmas=apkmas
562  END
563  FUNCTION lunpik(ID,ISGN)
564  COMMON / taukle / bra1,brk0,brk0b,brks
565  REAL*4 bra1,brk0,brk0b,brks
566  REAL*4 xio(1)
567  ident=id*isgn
568 #if defined (ALEPH)
569  IF (ident.EQ. 1) THEN
570  ipkdef= 211
571  ELSEIF (ident.EQ.-1) THEN
572  ipkdef=-211
573  ELSEIF (ident.EQ. 2) THEN
574  ipkdef= 111
575  ELSEIF (ident.EQ.-2) THEN
576  ipkdef= 111
577  ELSEIF (ident.EQ. 3) THEN
578  ipkdef= 321
579  ELSEIF (ident.EQ.-3) THEN
580  ipkdef=-321
581 #else
582  IF (ident.EQ. 1) THEN
583  ipkdef=-211
584  ELSEIF (ident.EQ.-1) THEN
585  ipkdef= 211
586  ELSEIF (ident.EQ. 2) THEN
587  ipkdef=111
588  ELSEIF (ident.EQ.-2) THEN
589  ipkdef=111
590  ELSEIF (ident.EQ. 3) THEN
591  ipkdef=-321
592  ELSEIF (ident.EQ.-3) THEN
593  ipkdef= 321
594 #endif
595  ELSEIF (ident.EQ. 4) THEN
596 *
597 * K0 --> K0_LONG (IS 130) / K0_SHORT (IS 310) = 1/1
598  CALL ranmar(xio,1)
599  IF (xio(1).GT.brk0) THEN
600  ipkdef= 130
601  ELSE
602  ipkdef= 310
603  ENDIF
604  ELSEIF (ident.EQ.-4) THEN
605 *
606 * K0B--> K0_LONG (IS 130) / K0_SHORT (IS 310) = 1/1
607  CALL ranmar(xio,1)
608  IF (xio(1).GT.brk0b) THEN
609  ipkdef= 130
610  ELSE
611  ipkdef= 310
612  ENDIF
613  ELSEIF (ident.EQ. 8) THEN
614  ipkdef= 22
615  ELSEIF (ident.EQ.-8) THEN
616  ipkdef= 22
617  ELSEIF (ident.EQ. 9) THEN
618  ipkdef= 221
619  ELSEIF (ident.EQ.-9) THEN
620  ipkdef= 221
621  ELSE
622  print *, 'STOP IN IPKDEF, WRONG IDENT=',ident
623  stop
624  ENDIF
625  lunpik=ipkdef
626  END
627 
628 
629 #if defined (CLEO)
630 
631  SUBROUTINE taurdf(KTO)
632 C THIS ROUTINE CAN BE CALLED BEFORE ANY TAU+ OR TAU- EVENT IS GENERATED
633 C IT CAN BE USED TO GENERATE TAU+ AND TAU- SAMPLES OF DIFFERENT
634 C CONTENTS
635  COMMON / taukle / bra1,brk0,brk0b,brks
636  REAL*4 bra1,brk0,brk0b,brks
637  COMMON / taubra / gamprt(30),jlist(30),nchan
638  IF (kto.EQ.1) THEN
639 C ==================
640 C AJWMOD: Set the BRs for (A1+ -> rho+ pi0) and (K*+ -> K0 pi+)
641  bra1 = pkorb(4,1)
642  brks = pkorb(4,3)
643  brk0 = pkorb(4,5)
644  brk0b = pkorb(4,6)
645  ELSE
646 C ====
647 C AJWMOD: Set the BRs for (A1+ -> rho+ pi0) and (K*+ -> K0 pi+)
648  bra1 = pkorb(4,2)
649  brks = pkorb(4,4)
650  brk0 = pkorb(4,5)
651  brk0b = pkorb(4,6)
652  ENDIF
653 C =====
654  END
655 #else
656 
657  SUBROUTINE taurdf(KTO)
658 * THIS ROUTINE CAN BE CALLED BEFORE ANY TAU+ OR TAU- EVENT IS GENERATED
659 * IT CAN BE USED TO GENERATE TAU+ AND TAU- SAMPLES OF DIFFERENT
660 * CONTENTS
661  COMMON / taukle / bra1,brk0,brk0b,brks
662  REAL*4 bra1,brk0,brk0b,brks
663  COMMON / taubra / gamprt(30),jlist(30),nchan
664  IF (kto.EQ.1) THEN
665 * ==================
666 * LIST OF BRANCHING RATIOS
667  nchan = 19
668  DO 1 i = 1,30
669  IF (i.LE.nchan) THEN
670  jlist(i) = i
671  IF(i.EQ. 1) gamprt(i) = .0000
672  IF(i.EQ. 2) gamprt(i) = .0000
673  IF(i.EQ. 3) gamprt(i) = .0000
674  IF(i.EQ. 4) gamprt(i) = .0000
675  IF(i.EQ. 5) gamprt(i) = .0000
676  IF(i.EQ. 6) gamprt(i) = .0000
677  IF(i.EQ. 7) gamprt(i) = .0000
678  IF(i.EQ. 8) gamprt(i) = 1.0000
679  IF(i.EQ. 9) gamprt(i) = 1.0000
680  IF(i.EQ.10) gamprt(i) = 1.0000
681  IF(i.EQ.11) gamprt(i) = 1.0000
682  IF(i.EQ.12) gamprt(i) = 1.0000
683  IF(i.EQ.13) gamprt(i) = 1.0000
684  IF(i.EQ.14) gamprt(i) = 1.0000
685  IF(i.EQ.15) gamprt(i) = 1.0000
686  IF(i.EQ.16) gamprt(i) = 1.0000
687  IF(i.EQ.17) gamprt(i) = 1.0000
688  IF(i.EQ.18) gamprt(i) = 1.0000
689  IF(i.EQ.19) gamprt(i) = 1.0000
690  ELSE
691  jlist(i) = 0
692  gamprt(i) = 0.
693  ENDIF
694  1 CONTINUE
695 * --- COEFFICIENTS TO FIX RATIO OF:
696 * --- A1 3CHARGED/ A1 1CHARGED 2 NEUTRALS MATRIX ELEMENTS (MASLESS LIM.)
697 * --- PROBABILITY OF K0 TO BE KS
698 * --- PROBABILITY OF K0B TO BE KS
699 * --- RATIO OF COEFFICIENTS FOR K*--> K0 PI-
700 * --- ALL COEFFICENTS SHOULD BE IN THE RANGE (0.0,1.0)
701 * --- THEY MEANING IS PROBABILITY OF THE FIRST CHOICE ONLY IF ONE
702 * --- NEGLECTS MASS-PHASE SPACE EFFECTS
703  bra1=0.5
704  brk0=0.5
705  brk0b=0.5
706  brks=0.6667
707  ELSE
708 * ====
709 * LIST OF BRANCHING RATIOS
710  nchan = 19
711  DO 2 i = 1,30
712  IF (i.LE.nchan) THEN
713  jlist(i) = i
714  IF(i.EQ. 1) gamprt(i) = .0000
715  IF(i.EQ. 2) gamprt(i) = .0000
716  IF(i.EQ. 3) gamprt(i) = .0000
717  IF(i.EQ. 4) gamprt(i) = .0000
718  IF(i.EQ. 5) gamprt(i) = .0000
719  IF(i.EQ. 6) gamprt(i) = .0000
720  IF(i.EQ. 7) gamprt(i) = .0000
721  IF(i.EQ. 8) gamprt(i) = 1.0000
722  IF(i.EQ. 9) gamprt(i) = 1.0000
723  IF(i.EQ.10) gamprt(i) = 1.0000
724  IF(i.EQ.11) gamprt(i) = 1.0000
725  IF(i.EQ.12) gamprt(i) = 1.0000
726  IF(i.EQ.13) gamprt(i) = 1.0000
727  IF(i.EQ.14) gamprt(i) = 1.0000
728  IF(i.EQ.15) gamprt(i) = 1.0000
729  IF(i.EQ.16) gamprt(i) = 1.0000
730  IF(i.EQ.17) gamprt(i) = 1.0000
731  IF(i.EQ.18) gamprt(i) = 1.0000
732  IF(i.EQ.19) gamprt(i) = 1.0000
733  ELSE
734  jlist(i) = 0
735  gamprt(i) = 0.
736  ENDIF
737  2 CONTINUE
738 * --- COEFFICIENTS TO FIX RATIO OF:
739 * --- A1 3CHARGED/ A1 1CHARGED 2 NEUTRALS MATRIX ELEMENTS (MASLESS LIM.)
740 * --- PROBABILITY OF K0 TO BE KS
741 * --- PROBABILITY OF K0B TO BE KS
742 * --- RATIO OF COEFFICIENTS FOR K*--> K0 PI-
743 * --- ALL COEFFICENTS SHOULD BE IN THE RANGE (0.0,1.0)
744 * --- THEY MEANING IS PROBABILITY OF THE FIRST CHOICE ONLY IF ONE
745 * --- NEGLECTS MASS-PHASE SPACE EFFECTS
746  bra1=0.5
747  brk0=0.5
748  brk0b=0.5
749  brks=0.6667
750  ENDIF
751 * =====
752  END
753 #endif
754 
755  SUBROUTINE iniphx(XK00)
756 * ----------------------------------------------------------------------
757 * INITIALISATION OF PARAMETERS
758 * USED IN QED and/or GSW ROUTINES
759 * ----------------------------------------------------------------------
760  COMMON / qedprm /alfinv,alfpi,xk0
761  REAL*8 alfinv,alfpi,xk0
762  REAL*8 pi8,xk00
763 *
764  pi8 = 4.d0*datan(1.d0)
765  alfinv = 137.03604d0
766  alfpi = 1d0/(alfinv*pi8)
767  xk0=xk00
768  END
769 
770  SUBROUTINE inimas
771 C ----------------------------------------------------------------------
772 C INITIALISATION OF MASSES
773 C
774 C called by : KORALZ
775 C ----------------------------------------------------------------------
776  COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
777  * ,ampiz,ampi,amro,gamro,ama1,gama1
778  * ,amk,amkz,amkst,gamkst
779 *
780  REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
781  * ,ampiz,ampi,amro,gamro,ama1,gama1
782  * ,amk,amkz,amkst,gamkst
783 C
784 C IN-COMING / OUT-GOING FERMION MASSES
785  amtau = 1.7842
786 C --- tau mass must be the same as in the host program, what-so-ever
787  amtau = 1.777
788  amnuta = 0.010
789  amel = 0.0005111
790  amnue = 0.0
791  ammu = 0.105659
792  amnumu = 0.0
793 *
794 * MASSES USED IN TAU DECAYS
795 #if defined (CePeCe)
796  ampiz = 0.134964
797  ampi = 0.139568
798  amro = 0.773
799  gamro = 0.145
800 *C GAMRO = 0.666
801  ama1 = 1.251
802  gama1 = 0.599
803  amk = 0.493667
804  amkz = 0.49772
805  amkst = 0.8921
806  gamkst = 0.0513
807 #elif defined (CLEO)
808  ampiz = 0.134964
809  ampi = 0.139568
810  amro = 0.773
811  gamro = 0.145
812 *C GAMRO = 0.666
813  ama1 = 1.251
814  gama1 = 0.599
815  amk = 0.493667
816  amkz = 0.49772
817  amkst = 0.8921
818  gamkst = 0.0513
819 C
820 C
821 C IN-COMING / OUT-GOING FERMION MASSES
822 !! AMNUTA = PKORB(1,2)
823 !! AMNUE = PKORB(1,4)
824 !! AMNUMU = PKORB(1,6)
825 C
826 C MASSES USED IN TAU DECAYS Cleo settings
827 !! AMPIZ = PKORB(1,7)
828 !! AMPI = PKORB(1,8)
829 !! AMRO = PKORB(1,9)
830 !! GAMRO = PKORB(2,9)
831  ama1 = 1.275 !! PKORB(1,10)
832  gama1 = 0.615 !! PKORB(2,10)
833 !! AMK = PKORB(1,11)
834 !! AMKZ = PKORB(1,12)
835 !! AMKST = PKORB(1,13)
836 !! GAMKST = PKORB(2,13)
837 C
838 #elif defined (ALEPH)
839  ampiz = 0.134964
840  ampi = 0.139568
841  amro = 0.7714
842  gamro = 0.153
843 cam AMRO = 0.773
844 cam GAMRO = 0.145
845  ama1 = 1.251! PMAS(LUCOMP(ia1),1) ! AMA1 = 1.251
846  gama1 = 0.599! PMAS(LUCOMP(ia1),2) ! GAMA1 = 0.599
847  print *,'INIMAS a1 mass= ',ama1,gama1
848  amk = 0.493667
849  amkz = 0.49772
850  amkst = 0.8921
851  gamkst = 0.0513
852 #else
853 #endif
854 
855  RETURN
856  END
857  subroutine bostdq(idir,vv,pp,q)
858 * *******************************
859 c Boost along arbitrary vector v (see eg. J.D. Jacson, Classical
860 c Electrodynamics).
861 c Four-vector pp is boosted from an actual frame to the rest frame
862 c of the four-vector v (for idir=1) or back (for idir=-1).
863 c q is a resulting four-vector.
864 c Note: v must be time-like, pp may be arbitrary.
865 c
866 c Written by: Wieslaw Placzek date: 22.07.1994
867 c Last update: 3/29/95 by: M.S.
868 c
869  implicit DOUBLE PRECISION (a-h,o-z)
870  parameter(nout=6)
871  DOUBLE PRECISION v(4),p(4),q(4),pp(4),vv(4)
872  save
873 !
874  do 1 i=1,4
875  v(i)=vv(i)
876  1 p(i)=pp(i)
877  amv=(v(4)**2-v(1)**2-v(2)**2-v(3)**2)
878  if (amv.le.0d0) then
879  write(6,*) 'bosstv: warning amv**2=',amv
880  endif
881  amv=sqrt(abs(amv))
882  if (idir.eq.-1) then
883  q(4)=( p(1)*v(1)+p(2)*v(2)+p(3)*v(3)+p(4)*v(4))/amv
884  wsp =(q(4)+p(4))/(v(4)+amv)
885  elseif (idir.eq.1) then
886  q(4)=(-p(1)*v(1)-p(2)*v(2)-p(3)*v(3)+p(4)*v(4))/amv
887  wsp =-(q(4)+p(4))/(v(4)+amv)
888  else
889  write(nout,*)' >>> boostv: wrong value of idir = ',idir
890  endif
891  q(1)=p(1)+wsp*v(1)
892  q(2)=p(2)+wsp*v(2)
893  q(3)=p(3)+wsp*v(3)
894  end
895 
896 
897 
898 
899 
900 
901 
902 
903 
904 #if defined (ALEPH)
905  FUNCTION dilogy(X)
906 C *****************
907  IMPLICIT REAL*8(a-h,o-z)
908 CERN C304 VERSION 29/07/71 DILOG 59 C
909  z=-1.64493406684822
910  IF(x .LT.-1.0) go to 1
911  IF(x .LE. 0.5) go to 2
912  IF(x .EQ. 1.0) go to 3
913  IF(x .LE. 2.0) go to 4
914  z=3.2898681336964
915  1 t=1.0/x
916  s=-0.5
917  z=z-0.5* log(abs(x))**2
918  go to 5
919  2 t=x
920  s=0.5
921  z=0.
922  go to 5
923  3 dilogy=1.64493406684822
924  RETURN
925  4 t=1.0-x
926  s=-0.5
927  z=1.64493406684822 - log(x)* log(abs(t))
928  5 y=2.66666666666666 *t+0.66666666666666
929  b= 0.00000 00000 00001
930  a=y*b +0.00000 00000 00004
931  b=y*a-b+0.00000 00000 00011
932  a=y*b-a+0.00000 00000 00037
933  b=y*a-b+0.00000 00000 00121
934  a=y*b-a+0.00000 00000 00398
935  b=y*a-b+0.00000 00000 01312
936  a=y*b-a+0.00000 00000 04342
937  b=y*a-b+0.00000 00000 14437
938  a=y*b-a+0.00000 00000 48274
939  b=y*a-b+0.00000 00001 62421
940  a=y*b-a+0.00000 00005 50291
941  b=y*a-b+0.00000 00018 79117
942  a=y*b-a+0.00000 00064 74338
943  b=y*a-b+0.00000 00225 36705
944  a=y*b-a+0.00000 00793 87055
945  b=y*a-b+0.00000 02835 75385
946  a=y*b-a+0.00000 10299 04264
947  b=y*a-b+0.00000 38163 29463
948  a=y*b-a+0.00001 44963 00557
949  b=y*a-b+0.00005 68178 22718
950  a=y*b-a+0.00023 20021 96094
951  b=y*a-b+0.00100 16274 96164
952  a=y*b-a+0.00468 63619 59447
953  b=y*a-b+0.02487 93229 24228
954  a=y*b-a+0.16607 30329 27855
955  a=y*a-b+1.93506 43008 6996
956  dilogy=s*t*(a-b)+z
957  RETURN
958 C=======================================================================
959 C===================END OF CPC PART ====================================
960 C=======================================================================
961  END
962 #endif