EPOS/epos167/epos-ids-165.f

   1 c-----------------------------------------------------------------------
   2       subroutine iclass(id,icl)
   3 c-----------------------------------------------------------------------
   4 c      determines hadron class
   5 c-----------------------------------------------------------------------
   6       ida=iabs(id)
   7       if(ida.eq.0)then
   8        icl=2
   9       elseif(ida.eq.130.or.ida.eq.230.or.ida.eq.20)then
  10        icl=3
  11       elseif(ida.eq.140.or.ida.eq.240.or.ida.eq.340.or.ida.eq.441)then
  12        icl=4
  13       elseif(ida.ge.100.and.ida.le.999)then
  14        icl=1
  15       elseif(ida.ge.1000.and.ida.le.9999)then
  16        icl=2
  17       else
  18        stop'iclass: id not known'
  19       endif
  20       end
  21
  22 c-----------------------------------------------------------------------
  23       subroutine idchrg(id,chrg)
  24 c     computes charge of particle with ident code id
  25 c     ichrg must be dimensioned nqlep+12
  26 c-----------------------------------------------------------------------
  27       dimension ichrg(53),ifl(3)
  28       data ichrg/0,2,-1,-1,2,-1,2,-1,2,0,0,0,-3,0,-3,0,-3,1,1,2,2*0
  29      *,2,-1,-1,2,-1,2,-1,2,0,0,0,-3,0,-3,0,-3,0,-3,3,0
  30      *,3,0,0,0,3,3,3,6,6,6,0/
  31       idabs=iabs(id)
  32       call idflav(id,ifl(1),ifl(2),ifl(3),jspin,ind)
  33       if(idabs.lt.100) goto 200
  34       isum=0
  35       do 100 i=1,3
  36       if(abs(ifl(i)).gt.52)goto 100
  37       isum=isum+ichrg(iabs(ifl(i))+1)*isign(1,ifl(i))
  38   100 continue
  39       chrg=isum/3.
  40       return
  41 200   chrg=ichrg(ind+1)*isign(1,id)
  42       chrg=chrg/3.
  43       return
  44       end
  45
  46 c-----------------------------------------------------------------------
  47       subroutine idcomk(ic)
  48 c     compactifies ic
  49 c-----------------------------------------------------------------------
  50       parameter (nflav=6)
  51       integer ic(2),icx(2),jc(nflav,2)
  52       call idcomp(ic,icx,jc,1)
  53       ic(1)=icx(1)
  54       ic(2)=icx(2)
  55       return
  56       end
  57
  58 cc-----------------------------------------------------------------------
  59 c      subroutine idcomi(ic,icx)
  60 cc     compactifies ic
  61 cc-----------------------------------------------------------------------
  62 c      parameter (nflav=6)
  63 c      integer ic(2),icx(2),jc(nflav,2)
  64 c      call idcomp(ic,icx,jc,1)
  65 c      return
  66 c      end
  67 c
  68 c-----------------------------------------------------------------------
  69       subroutine idcomj(jc)
  70 c     compactifies jc
  71 c-----------------------------------------------------------------------
  72       parameter (nflav=6)
  73       integer ic(2),icx(2),jc(nflav,2)
  74       call idcomp(ic,icx,jc,2)
  75       return
  76       end
  77
  78 c-----------------------------------------------------------------------
  79       subroutine idcomp(ic,icx,jc,im)
  80 c-----------------------------------------------------------------------
  81 c     compactifies ic,jc
  82 c     input: im (1 or 2)
  83 c            ic (if im=1)
  84 c            jc (if im=2)
  85 c     output: icx (if im=1)
  86 c             jc
  87 c-----------------------------------------------------------------------
  88       parameter (nflav=6)
  89       integer ic(2),icx(2),jc(nflav,2)
  90       if(im.eq.1)call iddeco(ic,jc)
  91       icx(1)=0
  92       icx(2)=0
  93            do n=1,nflav
  94            do j=1,2
  95       if(jc(n,j).ne.0)goto1
  96            enddo
  97            enddo
  98       return
  99 1     continue
 100       nq=0
 101       na=0
 102            do n=1,nflav
 103       nq=nq+jc(n,1)
 104       na=na+jc(n,2)
 105            enddo
 106       l=0
 107            do n=1,nflav
 108       k=min0(jc(n,1),jc(n,2))
 109       if(nq.eq.1.and.na.eq.1)k=0
 110       jc(n,1)=jc(n,1)-k
 111       jc(n,2)=jc(n,2)-k
 112       if(jc(n,1).lt.0.or.jc(n,2).lt.0)
 113      *call utstop('idcomp: jc negative&')
 114       l=l+jc(n,1)+jc(n,2)
 115            enddo
 116            if(l.eq.0)then
 117       jc(1,1)=1
 118       jc(1,2)=1
 119            endif
 120            if(im.eq.1)then
 121       call idenco(jc,icx,ireten)
 122       if(ireten.eq.1)call utstop('idcomp: idenco ret code = 1&')
 123            endif
 124       return
 125       end
 126
 127 c-----------------------------------------------------------------------
 128       subroutine iddeco(ic,jc)
 129 c     decode particle id
 130 c-----------------------------------------------------------------------
 131       parameter (nflav=6)
 132       integer jc(nflav,2),ic(2)
 133       ici=ic(1)
 134       jc(6,1)=mod(ici,10)
 135       jc(5,1)=mod(ici/10,10)
 136       jc(4,1)=mod(ici/100,10)
 137       jc(3,1)=mod(ici/1000,10)
 138       jc(2,1)=mod(ici/10000,10)
 139       jc(1,1)=mod(ici/100000,10)
 140       ici=ic(2)
 141       jc(6,2)=mod(ici,10)
 142       jc(5,2)=mod(ici/10,10)
 143       jc(4,2)=mod(ici/100,10)
 144       jc(3,2)=mod(ici/1000,10)
 145       jc(2,2)=mod(ici/10000,10)
 146       jc(1,2)=mod(ici/100000,10)
 147       return
 148       end
 149
 150 c-----------------------------------------------------------------------
 151       subroutine idenco(jc,ic,ireten)
 152 c     encode particle id
 153 c-----------------------------------------------------------------------
 154       parameter (nflav=6)
 155       integer jc(nflav,2),ic(2)
 156       ireten=0
 157       ic(1)=0
 158       do 20 i=1,nflav
 159       if(jc(i,1).ge.10)goto22
 160 20    ic(1)=ic(1)+jc(i,1)*10**(nflav-i)
 161       ic(2)=0
 162       do 21 i=1,nflav
 163       if(jc(i,2).ge.10)goto22
 164 21    ic(2)=ic(2)+jc(i,2)*10**(nflav-i)
 165       return
 166 22    ireten=1
 167       ic(1)=0
 168       ic(2)=0
 169       return
 170       end
 171
 172 c-----------------------------------------------------------------------
 173       subroutine idenct(jc,id,ib1,ib2,ib3,ib4)
 174 c     encode particle id
 175 c-----------------------------------------------------------------------
 176       parameter (nflav=6)
 177       integer jc(nflav,2),ic(2)
 178
 179       do 40 nf=1,nflav
 180       do 40 ij=1,2
 181       if(jc(nf,ij).ge.10)id=7*10**8
 182 40    continue
 183            if(id/10**8.ne.7)then
 184       call idenco(jc,ic,ireten)
 185       if(ireten.eq.1)call utstop('idenct: idenco ret code = 1&')
 186       if(mod(ic(1),100).ne.0.or.mod(ic(2),100).ne.0)then
 187       id=9*10**8
 188       else
 189       id=8*10**8+ic(1)*100+ic(2)/100
 190       endif
 191            else
 192       call idtrbi(jc,ib1,ib2,ib3,ib4)
 193       id=id
 194      *+mod(jc(1,1)+jc(2,1)+jc(3,1)+jc(4,1),10**4)*10**4
 195      *+mod(jc(1,2)+jc(2,2)+jc(3,2)+jc(4,2),10**4)
 196            endif
 197       return
 198       end
 199
 200 c-----------------------------------------------------------------------
 201       subroutine idflav(id,ifl1,ifl2,ifl3,jspin,index)
 202 c     unpacks the ident code id=+/-ijkl
 203 c
 204 c          mesons--
 205 c          i=0, j<=k, +/- is sign for j
 206 c          id=110 for pi0, id=220 for eta, etc.
 207 c
 208 c          baryons--
 209 c          i<=j<=k in general
 210 c          j<i<k for second state antisymmetric in (i,j), eg. l = 2130
 211 c
 212 c          other--
 213 c          id=1,...,6 for quarks
 214 c          id=9 for gluon
 215 c          id=10 for photon
 216 c          id=11,...,16 for leptons
 217 c          i=17 for deuteron
 218 c          i=18 for triton
 219 c          i=19 for alpha
 220 c          id=20 for ks, id=-20 for kl
 221 c
 222 c          i=21...26 for scalar quarks
 223 c          i=29 for gluino
 224 c
 225 c          i=30 for h-dibaryon
 226 c
 227 c          i=31...36 for scalar leptons
 228 c          i=39 for wino
 229 c          i=40 for zino
 230 c
 231 c          id=80 for w+
 232 c          id=81,...,83 for higgs mesons (h0, H0, A0, H+)
 233 c          id=84,...,87 for excited bosons (Z'0, Z''0, W'+)
 234 c          id=90 for z0
 235 c
 236 c          diquarks--
 237 c          id=+/-ij00, i<j for diquark composed of i,j.
 238 c
 239 c
 240 c          index is a sequence number used internally
 241 c          (index=-1 if id doesn't exist)
 242 c
 243 c-----------------------------------------------------------------------
 244       parameter ( nqlep=41,nmes=2)
 245       ifl1=0
 246       ifl2=0
 247       ifl3=0
 248       jspin=0
 249       idabs=iabs(id)
 250       i=idabs/1000
 251       j=mod(idabs/100,10)
 252       k=mod(idabs/10,10)
 253       jspin=mod(idabs,10)
 254       if(id.ge.10000) goto 400
 255       if(id.ne.0.and.mod(id,100).eq.0) goto 300
 256       if(j.eq.0) goto 200
 257       if(i.eq.0) goto 100
 258 c          baryons
 259 c          only x,y baryons are qqx, qqy, q=u,d,s.
 260       ifl1=isign(i,id)
 261       ifl2=isign(j,id)
 262       ifl3=isign(k,id)
 263       if(k.le.6) then
 264         index=max0(i-1,j-1)**2+i+max0(i-j,0)+(k-1)*k*(2*k-1)/6
 265      1  +109*jspin+36*nmes+nqlep+11
 266       else
 267         index=max0(i-1,j-1)**2+i+max0(i-j,0)+9*(k-7)+91
 268      1  +109*jspin+36*nmes+nqlep+11
 269       endif
 270       return
 271 c          mesons
 272 100   continue
 273       ifl1=0
 274       ifl2=isign(j,id)
 275       ifl3=isign(k,-id)
 276       index=j+k*(k-1)/2+36*jspin+nqlep
 277       index=index+11
 278       return
 279 c          quarks, leptons, etc
 280 200   continue
 281       ifl1=0
 282       ifl2=0
 283       ifl3=0
 284       jspin=0
 285       index=idabs
 286       if(idabs.lt.20) return
 287 c          define index=20 for ks, index=21 for kl
 288       index=idabs+1
 289       if(id.eq.20) index=20
 290 c          index=nqlep+1,...,nqlep+11 for w+, higgs, z0
 291       if(idabs.lt.80) return
 292       index=nqlep+idabs-79
 293       return
 294 300   ifl1=isign(i,id)
 295       ifl2=isign(j,id)
 296       ifl3=0
 297       jspin=0
 298       index=0
 299       return
 300  400  index=-1
 301       return
 302       end
 303
 304 c-----------------------------------------------------------------------
 305       subroutine idqufl(n,id,nqu,nqd,nqs)
 306 c     unpacks the ident code of particle (n) and give the number of
 307 c     quarks of each flavour(only u,d,s)
 308 c-----------------------------------------------------------------------
 309       include 'epos.inc'
 310       integer jc(nflav,2),ic(2)
 311
 312       nqu=0
 313       nqd=0
 314       nqs=0
 315       if(iabs(id).ge.7.and.iabs(id).lt.100.and.iabs(id).ne.20)return
 316       if(iabs(id)/10.eq.11.or.iabs(id)/10.eq.22)return
 317       if(iabs(id).eq.20)then
 318         if(iorptl(n).gt.0.and.idptl(iorptl(n)).gt.0)then
 319           nqd=1
 320           nqs=-1
 321         elseif(iorptl(n).gt.0)then
 322           nqd=-1
 323           nqs=1
 324         else
 325           if(ish.ge.4)write(ifch,*)'Cannot count the number of quark'
 326         endif
 327         return
 328       endif
 329       if(id.ne.0.and.mod(id,100).eq.0.and.id.le.10**8) goto 300
 330       if(id/10**8.ne.7)then
 331         call idtr4(id,ic)
 332         call iddeco(ic,jc)
 333       else
 334         call idtrb(ibptl(1,n),ibptl(2,n),ibptl(3,n),ibptl(4,n),jc)
 335       endif
 336       nqu=jc(1,1)-jc(1,2)
 337       nqd=jc(2,1)-jc(2,2)
 338       nqs=jc(3,1)-jc(3,2)
 339       return
 340  300  i=iabs(id)/1000
 341       j=mod(iabs(id)/100,10)
 342       ifl1=isign(i,id)
 343       ifl2=isign(j,id)
 344       if(iabs(ifl1).eq.1)nqu=isign(1,ifl1)
 345       if(iabs(ifl1).eq.2)nqd=isign(1,ifl1)
 346       if(iabs(ifl1).eq.3)nqs=isign(1,ifl1)
 347       if(iabs(ifl2).eq.1)nqu=nqu+isign(1,ifl2)
 348       if(iabs(ifl2).eq.2)nqd=nqd+isign(1,ifl2)
 349       if(iabs(ifl2).eq.3)nqs=nqs+isign(1,ifl2)
 350 c      write(ifch,*)'id',id,ifl1,ifl2,nqu,nqd,nqs
 351       return
 352       end
 353
 354 c-----------------------------------------------------------------------
 355       function idlabl(id)
 356 c     returns the character*8 label for the particle id
 357 c-----------------------------------------------------------------------
 358       parameter ( nqlep=41,nmes=2)
 359 c
 360       character*8 idlabl
 361       character*8 llep,lmes0,lmes1,lbar0,labar0,lbar1,labar1
 362       character*8 lqq,laqq
 363       dimension llep(104)
 364       dimension lmes0(64),lmes1(64)
 365       dimension lbar0(109),labar0(109),lbar1(109),labar1(109)
 366       dimension lqq(21),laqq(21)
 367 c          diquark labels
 368       data lqq/
 369      1'uu0. ','ud0. ','dd0. ','us0. ','ds0. ','ss0. ','uc0. ','dc0. ',
 370      2'sc0. ','cc0. ','ub0. ','db0. ','sb0. ','cb0. ','bb0. ','ut0. ',
 371      3'dt0. ','st0. ','ct0. ','bt0. ','tt0. '/
 372       data laqq/
 373      1'auu0.','aud0.','add0.','aus0.','ads0.','ass0.','auc0.','adc0.',
 374      2'asc0.','acc0.','aub0.','adb0.','asb0.','acb0.','abb0.','aut0.',
 375      3'adt0.','ast0.','act0.','abt0.','att0.'/
 376 c          quark and lepton labels
 377       data llep/
 378      *'     ','up   ','ub   ','dn   ','db   ','st   ','sb   ','ch   ',
 379      *'cb   ','bt   ','bb   ','tp   ','tb   ','y    ','yb   ','x    ',
 380      *'xb   ','gl   ','err  ','gm   ','err  ','nue  ','anue ','e-   ',
 381      *'e+   ','num  ','anum ','mu-  ','mu+  ','nut  ','anut ','tau- ',
 382      *'tau+ ','deut ','adeut','trit ','atrit','alph ','aalph','ks   ',
 383      *'err  ','err  ','kl   ',
 384      *'upss ','ubss ','dnss ','dbss ','stss ','sbss ','chss ','cbss ',
 385      *'btss ','bbss ','tpss ','tbss ','err  ','err  ','err  ','err  ',
 386      *'glss ','err  ','hdiba','err  ','ness ','aness','e-ss ','e+ss ',
 387      *'nmss ','anmss','mu-ss','mu+ss','ntss ','antss','t-ss ','t+ss ',
 388      *'err  ','err  ','err  ','err  ','w+ss ','w-ss ','z0ss ','err  ',
 389      *'w+   ','w-   ','h0   ','ah0  ','H0   ','aH0  ','A0   ','aA0  ',
 390      *'H+   ','H-   ','Zp0  ','aZp0 ','Zpp0 ','aZpp0','Wp+  ','Wp-  ',
 391      *'err  ','err  ','err  ','err  ','z0   '/
 392 c          0- meson labels
 393       data lmes0/
 394      1'pi0  ','pi+  ','eta  ','pi-  ','k+   ','k0   ','etap ','ak0  ',
 395      2'k-   ','ad0  ','d-   ','f-   ','etac ','f+   ','d+   ','d0   ',
 396      2'ub.  ','db.  ','sb.  ','cb.  ','bb.  ','bc.  ','bs.  ','bd.  ',
 397      3'bu.  ','ut.  ','dt.  ','st.  ','ct.  ','bt.  ','tt.  ','tb.  ',
 398      4'tc.  ','ts.  ','td.  ','tu.  ','uy.  ','dy.  ','sy.  ','cy.  ',
 399      5'by.  ','ty.  ','yy.  ','yt.  ','yb.  ','yc.  ','ys.  ','yd.  ',
 400      6'yu.  ','ux.  ','dx.  ','sx.  ','cx.  ','bx.  ','tx.  ','yx.  ',
 401      7'xx.  ','xy.  ','xt.  ','xb.  ','xc.  ','xs.  ','xd.  ','xu.  '/
 402 c          1- meson labels
 403       data lmes1/
 404      1'rho0 ','rho+ ','omeg ','rho- ','k*+  ','k*0  ','phi  ','ak*0 ',
 405      2'k*-  ','ad*0 ','d*-  ','f*-  ','jpsi ','f*+  ','d*+  ','d*0  ',
 406      3'ub*  ','db*  ','sb*  ','cb*  ','upsl ','bc*  ','bs*  ','bd*  ',
 407      4'bu*  ','ut*  ','dt*  ','st*  ','ct*  ','bt*  ','tt*  ','tb*  ',
 408      5'tc*  ','ts*  ','td*  ','tu*  ','uy*  ','dy*  ','sy*  ','cy*  ',
 409      6'by*  ','ty*  ','yy*  ','yt*  ','yb*  ','yc*  ','ys*  ','yd*  ',
 410      7'yu*  ','ux*  ','dx*  ','sx*  ','cx*  ','bx*  ','tx*  ','yx*  ',
 411      8'xx*  ','xy*  ','xt*  ','xb*  ','xc*  ','xs*  ','xd*  ','xu*  '/
 412 c          1/2+ baryon labels
 413       data lbar0/
 414      1'err  ','p    ','n    ','err  ','err  ','s+   ','s0   ','s-   ',
 415      2'l    ','xi0  ','xi-  ','err  ','err  ','err  ','sc++ ','sc+  ',
 416      3'sc0  ','lc+  ','usc. ','dsc. ','ssc. ','sdc. ','suc. ','ucc. ',
 417      4'dcc. ','scc. ','err  ','err  ','err  ','err  ','uub. ','udb. ',
 418      5'ddb. ','dub. ','usb. ','dsb. ','ssb. ','sdb. ','sub. ','ucb. ',
 419      6'dcb. ','scb. ','ccb. ','csb. ','cdb. ','cub. ','ubb. ','dbb. ',
 420      7'sbb. ','cbb. ','err  ','err  ','err  ','err  ','err  ','utt. ',
 421      8'udt. ','ddt. ','dut. ','ust. ','dst. ','sst. ','sdt. ','sut. ',
 422      9'uct. ','dct. ','sct. ','cct. ','cst. ','cdt. ','cut. ','ubt. ',
 423      1'dbt. ','sbt. ','cbt. ','bbt. ','bct. ','bst. ','bdt. ','but. ',
 424      2'utt. ','dtt. ','stt. ','ctt. ','btt. ','err  ','err  ','err  ',
 425      3'err  ','err  ','err  ','uuy. ','udy. ','ddy. ','duy. ','usy. ',
 426      4'dsy. ','ssy. ','sdy. ','suy. ','uux. ','udx. ','ddx. ','dux. ',
 427      5'usx. ','dsx. ','ssx. ','sdx. ','sux. '/
 428       data labar0/
 429      1'err  ','ap   ','an   ','err  ','err  ','as-  ','as0  ','as+  ',
 430      2'al   ','axi0 ','axi+ ','err  ','err  ','err  ','asc--','asc- ',
 431      3'asc0 ','alc- ','ausc.','adsc.','assc.','asdc.','asuc.','aucc.',
 432      4'adcc.','ascc.','err  ','err  ','err  ','err  ','auub.','audb.',
 433      5'addb.','adub.','ausb.','adsb.','assb.','asdb.','asub.','aucb.',
 434      6'adcb.','ascb.','accb.','acsb.','acdb.','acub.','aubb.','adbb.',
 435      7'asbb.','acbb.','err  ','err  ','err  ','err  ','err  ','autt.',
 436      8'audt.','addt.','adut.','aust.','adst.','asst.','asdt.','asut.',
 437      9'auct.','adct.','asct.','acct.','acst.','acdt.','acut.','aubt.',
 438      1'adbt.','asbt.','acbt.','abbt.','abct.','abst.','abdt.','abut.',
 439      2'autt.','adtt.','astt.','actt.','abtt.','err  ','err  ','err  ',
 440      3'err  ','err  ','err  ','auuy.','audy.','addy.','aduy.','ausy.',
 441      4'adsy.','assy.','asdy.','asuy.','auux.','audx.','addx.','adux.',
 442      5'ausx.','adsx.','assx.','asdx.','asux.'/
 443 c          3/2+ baryon labels
 444       data lbar1/
 445      1'dl++ ','dl+  ','dl0  ','dl-  ','err  ','s*+  ','s*0  ','s*-  ',
 446      2'err  ','xi*0 ','xi*- ','om-  ','err  ','err  ','uuc* ','udc* ',
 447      3'ddc* ','err  ','usc* ','dsc* ','ssc* ','err  ','err  ','ucc* ',
 448      4'dcc* ','scc* ','ccc* ','err  ','err  ','err  ','uub* ','udb* ',
 449      5'ddb* ','err  ','usb* ','dsb* ','ssb* ','err  ','err  ','ucb* ',
 450      6'dcb* ','scb* ','ccb* ','err  ','err  ','err  ','ubb* ','dbb* ',
 451      7'sbb* ','cbb* ','bbb* ','err  ','err  ','err  ','err  ','utt* ',
 452      8'udt* ','ddt* ','err  ','ust* ','dst* ','sst* ','err  ','err  ',
 453      9'uct* ','dct* ','sct* ','cct* ','err  ','err  ','err  ','ubt* ',
 454      1'dbt* ','sbt* ','cbt* ','bbt* ','err  ','err  ','err  ','err  ',
 455      2'utt* ','dtt* ','stt* ','ctt* ','btt* ','ttt* ','err  ','err  ',
 456      3'err  ','err  ','err  ','uuy* ','udy* ','ddy* ','err  ','usy* ',
 457      4'dsy* ','ssy* ','err  ','err  ','uux* ','udx* ','ddx* ','err  ',
 458      5'usx* ','dsx* ','ssx* ','err  ','err  '/
 459       data labar1/
 460      1'adl--','adl- ','adl0 ','adl+ ','err  ','as*- ','as*0 ','as*+ ',
 461      2'err  ','axi*0','axi*+','aom+ ','err  ','err  ','auuc*','audc*',
 462      3'addc*','err  ','ausc*','adsc*','assc*','err  ','err  ','aucc*',
 463      4'adcc*','ascc*','accc*','err  ','err  ','err  ','auub*','audb*',
 464      5'addb*','err  ','ausb*','adsb*','assb*','err  ','err  ','aucb*',
 465      6'adcb*','ascb*','accb*','err  ','err  ','err  ','aubb*','adbb*',
 466      7'asbb*','acbb*','abbb*','err  ','err  ','err  ','err  ','autt*',
 467      8'audt*','addt*','err  ','aust*','adst*','asst*','err  ','err  ',
 468      9'auct*','adct*','asct*','acct*','err  ','err  ','err  ','aubt*',
 469      1'adbt*','asbt*','acbt*','abbt*','err  ','err  ','err  ','err  ',
 470      2'autt*','adtt*','astt*','actt*','abtt*','attt*','err  ','err  ',
 471      3'err  ','err  ','err  ','auuy*','audy*','addy*','err  ','ausy*',
 472      4'adsy*','assy*','err  ','err  ','auux*','audx*','addx*','err  ',
 473      5'ausx*','adsx*','assx*','err  ','err  '/
 474 c          entry
 475       call idflav(id,ifl1,ifl2,ifl3,jspin,ind)
 476       if(iabs(id).lt.100) goto200
 477       if(iabs(id).lt.1000) goto100
 478       if(id.ne.0.and.mod(id,100).eq.0) goto300
 479 c          baryons
 480       ind=ind-109*jspin-36*nmes-nqlep
 481       ind=ind-11
 482       if(jspin.eq.0.and.id.gt.0) idlabl=lbar0(ind)
 483       if(jspin.eq.0.and.id.lt.0) idlabl=labar0(ind)
 484       if(jspin.eq.1.and.id.gt.0) idlabl=lbar1(ind)
 485       if(jspin.eq.1.and.id.lt.0) idlabl=labar1(ind)
 486       return
 487 c          mesons
 488 100   continue
 489       i=max0(ifl2,ifl3)
 490       j=-min0(ifl2,ifl3)
 491       ind=max0(i-1,j-1)**2+i+max0(i-j,0)
 492       if(jspin.eq.0) idlabl=lmes0(ind)
 493       if(jspin.eq.1) idlabl=lmes1(ind)
 494       return
 495 c          quarks, leptons, etc.
 496 200   continue
 497       ind=2*ind
 498       if(id.le.0) ind=ind+1
 499       idlabl=llep(ind)
 500       return
 501 300   i=iabs(ifl1)
 502       j=iabs(ifl2)
 503       ind=i+j*(j-1)/2
 504       if(id.gt.0) idlabl=lqq(ind)
 505       if(id.lt.0) idlabl=laqq(ind)
 506       return
 507       end
 508
 509 c-----------------------------------------------------------------------
 510       subroutine idmass(idi,amass)
 511 c     returns the mass of the particle with ident code id.
 512 c     (deuteron, triton and alpha mass come from Gheisha ???)
 513 c-----------------------------------------------------------------------
 514       dimension ammes0(15),ammes1(15),ambar0(30),ambar1(30)
 515       dimension amlep(52)
 516       parameter ( nqlep=41,nmes=2)
 517 c-c   data amlep/.3,.3,.5,1.6,4.9,30.,-1.,-1.,0.,0.,
 518       data amlep/.005,.009,.180,1.6,4.9,170.,-1.,-1.,0.,0.,0.,.511003e-3
 519      *     ,0.,.105661,0.,1.807,1.87656,2.8167,3.755,.49767,.49767,
 520      *     100.3,100.3,100.5,101.6,104.9,130.,2*-1.,100.,0.,
 521      *     100.,100.005,100.,100.1,100.,101.8,2*-1.,100.,100.,
 522      *     11*0./
 523 c          0- meson mass table
 524       data ammes0/.13496,.13957,.5488,.49367,.49767,.9576,1.8633
 525      1     ,1.8683,2.030,2.976,5.279,5.279,5.369,6.5940,9.460/
 526 c     1- meson mass table
 527       data ammes1/.770,.770,.7826,.8881,.8922,1.0196,2.006,2.0086
 528      1     ,2.140,3.097,5.325,5.325,5.507,6.602,9.859/
 529 c     1/2+ baryon mass table
 530       data ambar0/-1.,.93828,.93957,2*-1.,1.1894,1.1925,1.1974
 531      1     ,1.1156,1.3149,1.3213,3*-1.
 532      $     ,2.453               !15          sigma_c++!
 533      $     ,2.454               !            sigma_c+
 534      $     ,2.452               !            sigma_c0
 535      $     ,2.285               !            lambda_c+
 536      2     ,2.466               !19  1340   !Xi_c+
 537      $     ,2.50                !20  2340   !Xi_c0
 538      $     ,2.704               !21  3340   !
 539      $     ,2.5                 !22  3240
 540      $     ,2.5                 !23  3140
 541      $     ,3.55                !24  1440
 542      $     ,3.55                !25  2440
 543      $     ,3.70                !26  3440
 544      $     ,4*-1./
 545 c     3/2+ baryon mass table
 546       data ambar1/1.232,1.232,1.232,1.232,-1.,1.3823,1.3820
 547      1     ,1.3875,-1.,1.5318,1.5350,1.6722,2*-1.
 548      2     ,2.519               !15          sigma_c++
 549      $     ,2.52                !            sigma_c+
 550      $     ,2.517               !            sigma_c0
 551      $     ,-1.
 552      $     ,2.645
 553      $     ,2.644
 554      $     ,2.80
 555      $     ,2*-1.
 556      $     ,3.75
 557      $     ,3.75
 558      3     ,3.90
 559      $     ,4.80
 560      $     ,3*-1./
 561 c     entry
 562       id=idi
 563       amass=0.
 564 ctp060829      if(iabs(id).eq.30)then
 565 ctp060829        amass=amhdibar
 566 ctp060829        return
 567 ctp060829      endif
 568       if(idi.gt.10000)return
 569       if(idi.eq.0)id=1120                        !for air target
 570       call idflav(id,ifl1,ifl2,ifl3,jspin,ind)
 571       if(id.ne.0.and.mod(id,100).eq.0) goto400
 572       if(iabs(ifl1).ge.5.or.iabs(ifl2).ge.5.or.iabs(ifl3).ge.5)
 573      1     goto300
 574       if(ifl2.eq.0) goto200
 575       if(ifl1.eq.0) goto100
 576 c          baryons
 577       ind=ind-109*jspin-36*nmes-nqlep
 578       ind=ind-11
 579       amass=(1-jspin)*ambar0(ind)+jspin*ambar1(ind)
 580       return
 581 c          mesons
 582 100   continue
 583       ind=ind-36*jspin-nqlep
 584       ind=ind-11
 585       amass=(1-jspin)*ammes0(ind)+jspin*ammes1(ind)
 586       return
 587 c          quarks and leptons (+deuteron, triton, alpha, Ks and Kl)
 588 200   continue
 589       amass=amlep(ind)
 590       return
 591 c          b and t particles
 592 300   continue
 593       amass=amlep(iabs(ifl2))+amlep(iabs(ifl3))-.03+.04*jspin
 594       if(ifl1.ne.0) amass=amass+amlep(iabs(ifl1))
 595       return
 596 c          diquarks
 597 400   amass=amlep(iabs(ifl1))+amlep(iabs(ifl2))
 598       return
 599       end
 600
 601 cc-----------------------------------------------------------------------
 602 c      subroutine idmix(ic,jspin,icm,idm)
 603 cc     accounts for flavour mixing
 604 cc-----------------------------------------------------------------------
 605 c      parameter (nflav=6)
 606 c      real pmix1(3,2),pmix2(3,2)
 607 c      integer ic(2),icm(2)
 608 c      data pmix1/.25,.25,.5,0.,.5,1./,pmix2/.5,.5,1.,0.,0.,1./
 609 c      icm(1)=0
 610 c      icm(2)=0
 611 c      idm=0
 612 c      i=ic(1)
 613 c      if(i.ne.ic(2))return
 614 c      id=0
 615 c      if(i.eq.100000)id=1
 616 c      if(i.eq. 10000)id=2
 617 c      if(i.eq.  1000)id=3
 618 c      if(id.eq.0)return
 619 c      rnd=rangen()
 620 c      idm=int(pmix1(id,jspin+1)+rnd)+int(pmix2(id,jspin+1)+rnd)+1
 621 c      icm(1)=10**(nflav-idm)
 622 c      icm(2)=ic(1)
 623 c      idm=idm*100+idm*10+jspin
 624 c      return
 625 c      end
 626 c
 627 cc-----------------------------------------------------------------------
 628 c      subroutine idcleanjc(jc)
 629 cc-----------------------------------------------------------------------
 630 c      parameter (nflav=6)
 631 c      integer jc(nflav,2)
 632 c      ns=0
 633 c      do n=1,nflav
 634 c        jj=min(jc(n,1),jc(n,2))
 635 c        jc(n,1)=jc(n,1)-jj
 636 c        jc(n,2)=jc(n,2)-jj
 637 c        ns=ns+jc(n,1)+jc(n,2)
 638 c      enddo
 639 c      if(ns.eq.0)then
 640 c        jc(1,1)=1
 641 c        jc(1,2)=1
 642 c      endif
 643 c      end
 644 c
 645 c-----------------------------------------------------------------------
 646       subroutine idquacjc(jc,nqu,naq)
 647 c     returns quark content of jc
 648 c        jc(nflav,2) = jc-type particle identification code.
 649 c        nqu = # quarks
 650 c        naq = # antiquarks
 651 c-----------------------------------------------------------------------
 652       parameter (nflav=6)
 653       integer jc(nflav,2)
 654       nqu=0
 655       naq=0
 656       do 53 n=1,nflav
 657         nqu=nqu+jc(n,1)
 658 53      naq=naq+jc(n,2)
 659       return
 660       end
 661
 662 c-----------------------------------------------------------------------
 663       subroutine idquac(i,nq,ns,na,jc)
 664 c     returns quark content of ptl i from /cptl/ .
 665 c        nq = # quarks - # antiquarks
 666 c        ns = # strange quarks - # strange antiquarks
 667 c        na = # quarks + # antiquarks
 668 c        jc(nflav,2) = jc-type particle identification code.
 669 c-----------------------------------------------------------------------
 670       include 'epos.inc'
 671       integer jc(nflav,2),ic(2)
 672
 673       if(iabs(idptl(i)).eq.20)then
 674       idptl(i)=230
 675       if(rangen().lt..5)idptl(i)=-230
 676       goto9999
 677       endif
 678
 679       if(iabs(idptl(i)).lt.100)then
 680       nq=0
 681       ns=0
 682       do 1 n=1,nflav
 683       jc(n,1)=0
 684 1     jc(n,2)=0
 685       return
 686       endif
 687
 688 9999  if(idptl(i)/10**8.ne.7)then
 689       call idtr4(idptl(i),ic)
 690       call iddeco(ic,jc)
 691       else
 692       call idtrb(ibptl(1,i),ibptl(2,i),ibptl(3,i),ibptl(4,i),jc)
 693       endif
 694       na=0
 695       nq=0
 696       do 53 n=1,nflav
 697       na=na+jc(n,1)+jc(n,2)
 698 53    nq=nq+jc(n,1)-jc(n,2)
 699       ns=   jc(3,1)-jc(3,2)
 700       return
 701       end
 702
 703 cc-----------------------------------------------------------------------
 704 c      subroutine idquad(i,nq,na,jc)
 705 cc-----------------------------------------------------------------------
 706 cc     returns quark content of ptl i from /cptl/ .
 707 cc        nq = # quarks - # antiquarks
 708 cc        na = # quarks + # antiquarks
 709 cc        jc(nflav,2) = jc-type particle identification code.
 710 cc-----------------------------------------------------------------------
 711 c      include 'epos.inc'
 712 c      integer jc(nflav,2),ic(2)
 713 c
 714 c      id=idptl(i)
 715 c      if(iabs(id).eq.20)then
 716 c      id=230
 717 c      if(rangen().lt..5)id=-230
 718 c      goto9999
 719 c      endif
 720 c
 721 c      if(iabs(id).lt.100)then
 722 c      nq=0
 723 cc      ns=0
 724 c      do 1 n=1,nflav
 725 c      jc(n,1)=0
 726 c1     jc(n,2)=0
 727 c      return
 728 c      endif
 729 c
 730 c9999  if(id/10**8.ne.7)then
 731 c      call idtr4(id,ic)
 732 c      call iddeco(ic,jc)
 733 c      else
 734 c      call idtrb(ibptl(1,i),ibptl(2,i),ibptl(3,i),ibptl(4,i),jc)
 735 c      endif
 736 c      na=0
 737 c      nq=0
 738 c      do 53 n=1,nflav
 739 c      na=na+jc(n,1)+jc(n,2)
 740 c53    nq=nq+jc(n,1)-jc(n,2)
 741 cc      ns=   jc(3,1)-jc(3,2)
 742 c      return
 743 c      end
 744 c
 745 c-----------------------------------------------------------------------
 746       integer function idraflx(fc,icl,jc,j,c,iretso)
 747 c-----------------------------------------------------------------------
 748       include 'epos.inc'
 749       integer jc(nflav,2)
 750       character c*1
 751       rstrassave=rstras(icl)
 752       rstras(icl)=rstras(icl)*fc
 753       idraflx=idrafl(icl,jc,j,c,iretso)
 754       rstras(icl)=rstrassave
 755       end
 756
 757 c-----------------------------------------------------------------------
 758       integer function idrafl(icl,jc,j,c,iretso)
 759 c-----------------------------------------------------------------------
 760 c     returns random flavor,
 761 c     if : c='v' : according to jc
 762 c          c='s' : from sea
 763 c          c='d' : from sea for second quark in diquark
 764 c          c='c' : take out c quark first
 765 c     j=1 quark, j=2 antiquark, updates jc (if iretso=0)
 766 c     iretso=0  ok
 767 c           =1 : more than 9 quarks of same flavor attempted
 768 c-----------------------------------------------------------------------
 769       include 'epos.inc'
 770       integer jc(nflav,2)
 771       character c*1
 772
 773 c        write(ifch,*)'entry idrafl, j,c,jc: ',j,c,jc
 774
 775       if(c.eq.'s')then
 776         pud1=0.5
 777         pu=pud1**0.87                 !0.87=1./1.15) see ifl1 def in epos-fra
 778         pd=1.-pu
 779         ps=0.5
 780         pu=rstrau(icl)*pu
 781         pd=rstrad(icl)*pd
 782         ps=rstras(icl)*ps
 783         pc=0.
 784       elseif(c.eq.'d')then
 785         pud2=0.5
 786         pu=pud2**0.87   !0.87=1./1.15) see ifl1 def in epos-fra
 787         pd=1.-pu
 788         ps=puds*0.5
 789         pu=pu*rstrau(icl)
 790         pd=pd*rstrad(icl)
 791         ps=ps*rstras(icl)
 792         pc=0.
 793       elseif(c.eq.'v')then
 794         pu=float(jc(1,j))
 795         pd=float(jc(2,j))
 796         ps=float(jc(3,j))
 797         pc=float(jc(4,j))
 798       elseif(c.eq.'c')then
 799         pu=0.
 800         pd=0.
 801         ps=0.
 802         pc=1.
 803       else
 804         stop'idrafl: dunnowhatodo'
 805       endif
 806
 807 c      write(ifch,*)'idrafl',c,pu,pd,ps
 808
 809       s=pu+pd+ps+pc
 810       if(s.gt.0.)then
 811        r=rangen()*s
 812        if(r.gt.(pu+pd+ps).and.pc.gt.0d0)then
 813         i=4
 814        elseif(r.gt.(pu+pd).and.ps.gt.0d0)then
 815         i=3
 816        elseif(r.gt.pu.and.pd.gt.0d0)then
 817         i=2
 818        else
 819         i=1
 820        endif
 821       else
 822         i=1+int((2.+rstras(icl))*rangen())
 823       endif
 824       idrafl=i
 825
 826 c      write(ifch,*)'jc before updating',jc
 827 c      write(ifch,*)'i,j,jc',i,j,jc
 828
 829       call idsufl(i,j,jc,iretso)
 830
 831       if(iretso.ne.0.and.ish.ge.2)then
 832         call utmsg('idrafl')
 833         write(ifmt,*)'iret none 0 in idrafl',iretso
 834         write(ifch,*)'iret none 0 in idrafl',iretso
 835         call utmsgf
 836       endif
 837
 838 c      write(ifch,*)'jc after updating',jc
 839
 840       return
 841       end
 842
 843 c-----------------------------------------------------------------------
 844       integer function idraflz(jc,j)
 845 c-----------------------------------------------------------------------
 846       include 'epos.inc'
 847       integer jc(nflav,2)
 848
 849       pu=float(jc(1,j))
 850       pd=float(jc(2,j))
 851       ps=float(jc(3,j))
 852
 853       s=pu+pd+ps
 854       if(s.gt.0.)then
 855        r=rangen()*s
 856        if(r.gt.(pu+pd).and.ps.gt.0d0)then
 857         i=3
 858        elseif(r.gt.pu.and.pd.gt.0d0)then
 859         i=2
 860        else
 861         i=1
 862        endif
 863       else
 864        stop'in idraflz (1)                      '
 865       endif
 866       idraflz=i
 867
 868       if(jc(i,j).lt.1)stop'in idraflz (2)              '
 869       jc(i,j)=jc(i,j)-1
 870
 871       end
 872
 873 c-----------------------------------------------------------------------
 874       subroutine idsufl(i,j,jc,iretso)
 875 c-----------------------------------------------------------------------
 876 c subtract flavor i, j=1 quark, j=2 antiquark
 877 c add antiflavor if jc(i,j)=0
 878 c iretso=0  ok
 879 c       =1 : more than 9 quarks of same flavor attempted
 880 c-----------------------------------------------------------------------
 881       integer jc(6,2),ic(2)
 882
 883       if(jc(i,j).gt.0)then
 884        jc(i,j)=jc(i,j)-1
 885        call idenco(jc,ic,iret)
 886        if(ic(1).eq.0.and.ic(2).eq.0)then
 887          jc(i,j)=jc(i,j)+1
 888          if(jc(i,3-j).lt.9.and.iret.eq.0)then
 889            jc(i,3-j)=jc(i,3-j)+1
 890          else
 891            iretso=1
 892          endif
 893        endif
 894       else
 895         if(j.eq.1)then
 896           if(jc(i,2).lt.9)then
 897             jc(i,2)=jc(i,2)+1
 898           else
 899             iretso=1
 900           endif
 901         else
 902           if(jc(i,1).lt.9)then
 903             jc(i,1)=jc(i,1)+1
 904           else
 905             iretso=1
 906           endif
 907         endif
 908       endif
 909
 910       return
 911       end
 912
 913 c-----------------------------------------------------------------------
 914       subroutine idsufl2(i,j,jc,iret)
 915 c-----------------------------------------------------------------------
 916 c substract flavor i, by adding antiquark i, j=1 quark, j=2 antiquark
 917 c Can replace idsufl if we don't want to cancel quarks and antiquarks
 918 c-----------------------------------------------------------------------
 919       parameter(nflav=6)
 920       integer jc(nflav,2)
 921
 922       if(jc(i,3-j).lt.9)then
 923         jc(i,3-j)=jc(i,3-j)+1
 924       else
 925         iret=1
 926       endif
 927
 928       return
 929       end
 930
 931 cc-----------------------------------------------------------------------
 932 c      subroutine idchfl(jc1,jc2,iret)
 933 cc-----------------------------------------------------------------------
 934 cc checks whether jc1 and jc2 have the same number of quarks and antiquarks
 935 cc if yes: iret=0, if no: iret=1
 936 cc-----------------------------------------------------------------------
 937 c      integer jc1(6,2),jc2(6,2)
 938 c
 939 c      iret=0
 940 c
 941 c      do j=1,2
 942 c       n1=0
 943 c       n2=0
 944 c       do i=1,6
 945 c        n1=n1+jc1(i,j)
 946 c        n2=n2+jc2(i,j)
 947 c       enddo
 948 c       if(n1.ne.n2)then
 949 c        iret=1
 950 c        return
 951 c       endif
 952 c      enddo
 953 c
 954 c      end
 955 c
 956 c
 957 c-----------------------------------------------------------------------
 958       subroutine idres(idi,am,idr,iadj,iii)
 959 c     returns resonance id idr corresponding to mass am.
 960 c     performs mass adjustment, if necessary (if so iadj=1, 0 else)
 961 c     (only for mesons and baryons, error (stop) otherwise)
 962 c-----------------------------------------------------------------------
 963       include 'epos.inc'
 964       parameter (mxindx=1000,mxre=100,mxma=11,mxmx=6)
 965       common/crema/indx(mxindx),rema(mxre,mxma),rewi(mxre,mxma)
 966      *,idmx(mxma,mxmx),icre1(mxre,mxma),icre2(mxre,mxma)
 967       character cad*10
 968
 969       write(cad,'(i10)')idi
 970       iadj=0
 971       idr=0
 972       if(idi.eq.10)return
 973       if(abs(am).lt.1.e-5)am=1e-5
 974       id=idi
 975       ami=am
 976       if(am.lt.0.)then
 977         call idmass(id,am)
 978         iadj=1
 979         if(am.le.0.)then
 980         write(ifch,*)'*****  warning in idres (0): '
 981      *,'neg mass returned from idmass'
 982         write(ifch,*)'id,am(input):',idi,ami
 983         am=1e-5
 984         endif
 985       endif
 986
 987       if(id.eq.0)goto 9999
 988       if(abs(id).eq.20)id=sign(230,idi)
 989       m1=1
 990       if(iabs(id).ge.1000)m1=3
 991       m2=2
 992       if(iabs(id).ge.1000)m2=mxmx
 993       do 3 k=m1,m2
 994       do 3 m=2,mxma
 995         if(iabs(id).eq.idmx(m,k)) then
 996           id=idmx(1,k)*10*id/iabs(id)
 997           goto 43
 998         endif
 999  3    continue
1000  43   continue
1001       ix=iabs(id)/10
1002       if(ix.lt.1.or.ix.gt.mxindx)then
1003         call utstop('idres: ix out of range. id='//cad//'&')
1004       endif
1005       i=indx(ix)
1006       if(i.lt.1.or.i.gt.mxre)then
1007         write(ifch,*)'idres problem',id,am
1008         call utstop('idres: particle not in table&')
1009       endif
1010       do 1 j=1,mxma-1
1011       if(am.ge.rema(i,j).and.am.le.rema(i,j+1))then
1012       if(j-1.gt.9)call utstop('idres: spin > 9&')
1013       idr=id/10*10+(j-1)*id/iabs(id)
1014       goto 2
1015       endif
1016 1     continue
1017       goto 9999
1018 2     continue
1019
1020       do 4 k=1,mxmx
1021       if(ix.eq.idmx(1,k))then
1022       if(j.lt.1.or.j.gt.mxma-1)
1023      *call utstop('idres: index j out of range&')
1024       if(idmx(j+1,k).ne.0)idr=idmx(j+1,k)*id/iabs(id)
1025       endif
1026 4     continue
1027
1028       iy=mod(iabs(idr),10)
1029       if(iy.gt.maxres)then
1030       iadj=0
1031       idr=0
1032       goto 9999
1033       endif
1034
1035       if(iy.ne.0.and.iy.ne.1)goto 9999
1036
1037       call idmass(idr,am)
1038       if(am.lt.0.)then
1039       write(ifch,*)'*****  error in idres: '
1040      *,'neg mass returned from idmass'
1041       write(ifch,*)'id,am(input):',idi,ami
1042       write(ifch,*)'idr,am:',idr,am
1043       call utstop('idres: neg mass returned from idmass&')
1044       endif
1045       del=max(1.e-3,2.*rewi(i,j))
1046       if(abs(ami-am).gt.del)iadj=1
1047 c      write(ifch,*)'res:',id,idr,ami,am,rewi(i,j),iadj
1048
1049 9999  continue
1050       if(iii.eq.1)then
1051         if(idi.eq.221)stop'\n\n     STOP in idres (1) \n\n'
1052         if(idr.eq.221)stop'\n\n     STOP in idres (2) \n\n'
1053         if(idr.eq.111)then
1054           if(rangen().le.0.5)idr=221
1055           call idmass(idr,am)
1056         endif
1057       endif
1058       if(.not.(ish.ge.8))return
1059       write(ifch,*)'return from idres. id,ami,am,idr,iadj:'
1060       write(ifch,*)idi,ami,am,idr,iadj
1061       return
1062       end
1063
1064 c-----------------------------------------------------------------------
1065       subroutine idresi
1066 c-----------------------------------------------------------------------
1067 c  initializes /crema/
1068 c  width for 151, 251, 351 arbitrary (no data found) !!!!!!!!!!!
1069 c-----------------------------------------------------------------------
1070
1071       parameter (mxindx=1000,mxre=100,mxma=11,mxmx=6)
1072       common/crema/indx(mxindx),rema(mxre,mxma),rewi(mxre,mxma)
1073      *,idmx(mxma,mxmx),icre1(mxre,mxma),icre2(mxre,mxma)
1074       parameter (n=33)
1075       dimension remai(n,mxma),rewii(n,mxma),idmxi(mxma,mxmx)
1076      *,icrei(n,2*mxma)
1077
1078       data (idmxi(j,1),j=1,mxma)/ 11, 110, 111,   0,   0,   0,   0, 4*0/
1079       data (idmxi(j,2),j=1,mxma)/ 22, 220, 330, 331,   0,   0,   0, 4*0/
1080       data (idmxi(j,3),j=1,mxma)/123,2130,1230,1231,1233,1234,1235, 4*0/
1081       data (idmxi(j,4),j=1,mxma)/124,2140,1240,1241,   0,   0,   0, 4*0/
1082       data (idmxi(j,5),j=1,mxma)/134,3140,1340,1341,   0,   0,   0, 4*0/
1083       data (idmxi(j,6),j=1,mxma)/234,3240,2340,2341,   0,   0,   0, 4*0/
1084
1085       data ((icrei(k,m),m=1,2*mxma),k=1,10)/
1086      *111,000000, 9*300000,    11*0,
1087      *222,000000, 9*030000,    11*0,
1088      *112,       10*210000,    11*0,
1089      *122,       10*120000,    11*0,
1090      *113,       10*201000,    11*0,
1091      *223,       10*021000,    11*0,
1092      *123,       10*111000,    11*0,
1093      *133,       10*102000,    11*0,
1094      *233,       10*012000,    11*0,
1095      *333,000000, 9*003000,    11*0/
1096       data ((icrei(k,m),m=1,2*mxma),k=11,20)/
1097      *114,       10*200100,    11*0,
1098      *124,       10*110100,    11*0,
1099      *224,       10*020100,    11*0,
1100      *134,       10*101100,    11*0,
1101      *234,       10*011100,    11*0,
1102      *334,       10*002100,    11*0,
1103      *144,       10*100200,    11*0,
1104      *244,       10*010200,    11*0,
1105      *344,       10*001200,    11*0,
1106      *444,000000, 9*000300,    11*0/
1107       data ((icrei(k,m),m=1,2*mxma),k=21,29)/
1108      * 11,  10*100000,    0,   10*100000,
1109      * 22,  10*001000,    0,   10*001000,
1110      * 12,  10*100000,    0,   10*010000,
1111      * 13,  10*100000,    0,   10*001000,
1112      * 23,  10*010000,    0,   10*001000,
1113      * 14,  10*100000,    0,   10*000100,
1114      * 24,  10*010000,    0,   10*000100,
1115      * 34,  10*001000,    0,   10*000100,
1116      * 44,  10*000100,    0,   10*000100/
1117       data ((icrei(k,m),m=1,2*mxma),k=30,33)/
1118      * 15,  10*100000,    0,   10*000010,
1119      * 25,  10*010000,    0,   10*000010,
1120      * 35,  10*001000,    0,   10*000010,
1121      *  3,  10*222000,    0,   10*000010/
1122
1123       data ((remai(k,m),m=1,mxma),k=1,10)/
1124      *111.,0.000,1.425,1.660,1.825,2.000,0.000,0.000,0.000,0.000,0.000,
1125      *222.,0.000,1.425,1.660,1.825,2.000,0.000,0.000,0.000,0.000,0.000,
1126      *112.,1.080,1.315,1.485,1.575,1.645,1.685,1.705,1.825,2.000,0.000,
1127      *122.,1.080,1.315,1.485,1.575,1.645,1.685,1.705,1.825,2.000,0.000,
1128      *113.,1.300,1.500,1.700,1.850,2.000,0.000,0.000,0.000,0.000,0.000,
1129      *223.,1.300,1.500,1.700,1.850,2.000,0.000,0.000,0.000,0.000,0.000,
1130      *123.,1.117,1.300,1.395,1.465,1.540,1.655,1.710,1.800,1.885,2.000,
1131 c     *123.,1.154,1.288,1.395,1.463,1.560,1.630,1.710,1.800,1.885,2.000,
1132      *133.,1.423,2.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1133      *233.,1.428,2.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1134 c     *133.,1.423,1.638,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1135 c     *233.,1.427,1.634,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1136      *333.,0.000,2.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000/
1137       data ((remai(k,m),m=1,mxma),k=11,20)/
1138      *114.,2.530,2.730,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1139      *124.,2.345,2.530,2.730,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1140      *224.,2.530,2.730,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1141      *134.,2.450,2.600,2.800,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1142      *234.,2.450,2.600,2.800,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1143      *334.,2.700,2.900,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1144      *144.,3.650,3.850,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1145      *244.,3.650,3.850,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1146      *344.,3.800,4.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1147      *444.,0.000,5.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000/
1148       data ((remai(k,m),m=1,mxma),k=21,29)/
1149      * 11.,0.450,0.950,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1150      * 22.,0.750,0.965,1.500,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1151      * 12.,0.450,0.950,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1152      * 13.,0.500,1.075,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1153      * 23.,0.500,1.075,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1154      * 14.,1.935,2.150,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1155      * 24.,1.938,2.150,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1156      * 34.,2.085,2.370,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1157      * 44.,3.037,3.158,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000/
1158       data ((remai(k,m),m=1,mxma),k=30,33)/
1159      * 15.,5.302,5.348,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1160      * 25.,5.302,5.348,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1161      * 35.,5.390,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,
1162      *  3.,2.230,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000/
1163
1164       data ((rewii(k,m),m=1,mxma),k=1,5)/
1165      *111.,0.000e+00,0.115e+00,0.140e+00,0.250e+00,0.250e+00,
1166      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1167      *222.,0.000e+00,0.115e+00,0.140e+00,0.250e+00,0.250e+00,
1168      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1169      *112.,0.000e+00,0.115e+00,0.200e+00,0.140e+00,0.140e+00,
1170      *     0.145e+00,0.250e+00,0.140e+00,0.250e+00,0.000e+00,
1171      *122.,0.000e+00,0.115e+00,0.200e+00,0.140e+00,0.140e+00,
1172      *     0.145e+00,0.250e+00,0.140e+00,0.250e+00,0.000e+00,
1173      *113.,0.824e-14,0.036e+00,0.080e+00,0.100e+00,0.170e+00,
1174      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00/
1175       data ((rewii(k,m),m=1,mxma),k=6,10)/
1176      *223.,0.445e-14,0.039e+00,0.080e+00,0.100e+00,0.170e+00,
1177      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1178      *123.,0.250e-14,0.890e-05,0.036e+00,0.040e+00,0.016e+00,
1179      *     0.090e+00,0.080e+00,0.100e+00,0.145e+00,0.170e+00,
1180      *133.,0.227e-14,0.009e+00,0.000e+00,0.000e+00,0.000e+00,
1181      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1182      *233.,0.400e-14,0.010e+00,0.000e+00,0.000e+00,0.000e+00,
1183      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1184      *333.,0.000e+00,0.800e-14,0.000e+00,0.000e+00,0.000e+00,
1185      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00/
1186       data ((rewii(k,m),m=1,mxma),k=11,15)/
1187      *114.,0.400e-11,0.010e+00,0.000e+00,0.000e+00,0.000e+00,
1188      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1189      *124.,0.400e-11,0.400e-11,0.010e+00,0.000e+00,0.000e+00,
1190      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1191      *224.,0.400e-11,0.010e+00,0.010e+00,0.000e+00,0.000e+00,
1192      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1193      *134.,0.150e-11,0.400e-11,0.010e+00,0.000e+00,0.000e+00,
1194      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1195      *234.,0.150e-11,0.400e-11,0.010e+00,0.000e+00,0.000e+00,
1196      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00/
1197       data ((rewii(k,m),m=1,mxma),k=16,20)/
1198      *334.,0.400e-11,0.010e+00,0.010e+00,0.000e+00,0.000e+00,
1199      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1200      *144.,0.400e-11,0.010e+00,0.010e+00,0.000e+00,0.000e+00,
1201      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1202      *244.,0.400e-11,0.010e+00,0.010e+00,0.000e+00,0.000e+00,
1203      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1204      *344.,0.400e-11,0.010e+00,0.010e+00,0.000e+00,0.000e+00,
1205      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1206      *444.,0.400e-11,0.010e+00,0.010e+00,0.000e+00,0.000e+00,
1207      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00/
1208       data ((rewii(k,m),m=1,mxma),k=21,25)/
1209      * 11.,0.757e-08,0.153e+00,0.057e+00,0.000e+00,0.000e+00,
1210      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1211      * 22.,0.105e-05,0.210e-03,0.034e+00,0.004e+00,0.000e+00,
1212      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1213      * 12.,0.000e+00,0.153e+00,0.057e+00,0.000e+00,0.000e+00,
1214      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1215      * 13.,0.000e+00,0.051e+00,0.000e+00,0.000e+00,0.000e+00,
1216      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1217      * 23.,0.197e-02,0.051e+00,0.000e+00,0.000e+00,0.000e+00,
1218      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00/
1219       data ((rewii(k,m),m=1,mxma),k=26,29)/
1220      * 14.,0.154e-11,0.002e+00,0.000e+00,0.000e+00,0.000e+00,
1221      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1222      * 24.,0.615e-12,0.002e+00,0.000e+00,0.000e+00,0.000e+00,
1223      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1224      * 34.,0.150e-11,0.020e+00,0.000e+00,0.000e+00,0.000e+00,
1225      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1226      * 44.,0.010e+00,0.068e-03,0.000e+00,0.000e+00,0.000e+00,
1227      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00/
1228       data ((rewii(k,m),m=1,mxma),k=30,33)/
1229      * 15.,0.426e-12,0.010e+00,0.000e+00,0.000e+00,0.000e+00,
1230      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1231      * 25.,0.426e-12,0.010e+00,0.000e+00,0.000e+00,0.000e+00,
1232      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1233      * 35.,0.408e-12,0.010e+00,0.000e+00,0.000e+00,0.000e+00,
1234      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1235      *  3.,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,
1236      *     0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00/
1237
1238       do 3 i=1,mxindx
1239 3     indx(i)=0
1240       do 4 k=1,mxre
1241       do 4 m=1,mxma
1242 4     rema(k,m)=0
1243
1244       do 2 j=1,mxma
1245       do 2 i=1,mxmx
1246 2     idmx(j,i)=idmxi(j,i)
1247
1248       ntec=n
1249       if(ntec.gt.mxre)call utstop('idresi: dimension mxre too small&')
1250       do 1 k=1,n
1251       ix=nint(remai(k,1))
1252       ix2=nint(rewii(k,1))
1253       ix3=icrei(k,1)
1254       if(ix.ne.ix2)call utstop('idresi: ix /= ix2&')
1255       if(ix.ne.ix3)call utstop('idresi: ix /= ix3&')
1256       if(ix.lt.1.or.ix.gt.mxindx)
1257      *call utstop('idresi: ix out of range.&')
1258       indx(ix)=k
1259       rema(k,1)=0.
1260       rewi(k,1)=0.
1261       icre1(k,1)=0
1262       icre2(k,1)=0
1263       do 1 m=2,mxma
1264       rema(k,m)=remai(k,m)
1265       rewi(k,m)=rewii(k,m)
1266       icre1(k,m)=icrei(k,m)
1267 1     icre2(k,m)=icrei(k,mxma+m)
1268
1269       indx(33) =indx(22)
1270       indx(213)=indx(123)
1271       indx(214)=indx(124)
1272       indx(314)=indx(134)
1273       indx(324)=indx(234)
1274
1275       return
1276       end
1277
1278 cc-----------------------------------------------------------------------
1279 c      integer function idsgl(ic,gen,cmp)
1280 cc     returns 1 for singlets (qqq or qqbar) 0 else.
1281 cc-----------------------------------------------------------------------
1282 c      parameter (nflav=6)
1283 c      integer ic(2),jcx(nflav,2),icx(2)
1284 c      character gen*6,cmp*6
1285 c      idsgl=0
1286 c      if(cmp.eq.'cmp-ys')then
1287 c      call idcomi(ic,icx)
1288 c      else
1289 c      icx(1)=ic(1)
1290 c      icx(2)=ic(2)
1291 c      endif
1292 c      call iddeco(icx,jcx)
1293 c      nq=0
1294 c      na=0
1295 c      do 1 i=1,nflav
1296 c      nq=nq+jcx(i,1)
1297 c1     na=na+jcx(i,2)
1298 c      if(nq.eq.0.and.na.eq.0)return
1299 c      if(gen.eq.'gen-no')then
1300 c      if(nq.eq.3.and.na.eq.0.or.nq.eq.1.and.na.eq.1
1301 c     *.or.nq.eq.0.and.na.eq.3)idsgl=1
1302 c      elseif(gen.eq.'gen-ys')then
1303 c      if(mod(nq-na,3).eq.0)idsgl=1
1304 c      endif
1305 c      return
1306 c      end
1307 c
1308
1309 c-----------------------------------------------------------------------
1310       subroutine idtau(id,p4,p5,taugm)
1311 c     returns lifetime*gamma for id with energy p4, mass p5
1312 c-----------------------------------------------------------------------
1313       include 'epos.inc'
1314       parameter (mxindx=1000,mxre=100,mxma=11,mxmx=6)
1315       common/crema/indx(mxindx),rema(mxre,mxma),rewi(mxre,mxma)
1316      *,idmx(mxma,mxmx),icre1(mxre,mxma),icre2(mxre,mxma)
1317            if(iabs(id).lt.100.and.id.ne.20)then
1318       wi=0
1319            elseif(id.eq.20)then
1320       wi=.197/2.675e13
1321            elseif(id.eq.221)then
1322       wi=.197/0.008
1323            elseif(iabs(id).lt.1e8)then
1324       ix=iabs(id)/10
1325       if(ix.lt.1.or.ix.gt.mxindx)
1326      *call utstop('idtau: ix out of range.&')
1327       ii=indx(ix)
1328       jj=mod(iabs(id),10)+2
1329
1330       m1=1
1331       if(iabs(id).ge.1000)m1=3
1332       m2=2
1333       if(iabs(id).ge.1000)m2=mxmx
1334       do 75 imx=m1,m2
1335       do 75 ima=2,mxma
1336       if(iabs(id).eq.idmx(ima,imx))then
1337         jj=ima
1338         goto 75
1339       endif
1340 75    continue
1341  76   if(ii.lt.1.or.ii.gt.mxre.or.jj.lt.1.or.jj.gt.mxma)then
1342       write(ifch,*)'id,ii,jj:',id,'   ',ii,jj
1343       call utstop('idtau: ii or jj out of range&')
1344       endif
1345       wi=rewi(ii,jj)
1346            else
1347       tauz=taunll
1348 c-c   tauz=amin1(9./p5**2,tauz)
1349 c-c   tauz=amax1(.2,tauz)
1350       wi=.197/tauz
1351            endif
1352       if(wi.eq.0.)then
1353       tau=ainfin
1354       else
1355       tau=.197/wi
1356       endif
1357       if(p5.ne.0.)then
1358       gm=p4/p5
1359       else
1360       gm=ainfin
1361       endif
1362       if(tau.ge.ainfin.or.gm.ge.ainfin)then
1363       taugm=ainfin
1364       else
1365       taugm=tau*gm
1366       endif
1367       return
1368       end
1369
1370 c-----------------------------------------------------------------------
1371       subroutine idtr4(id,ic)
1372 c     transforms generalized paige_id -> werner_id  (for < 4 flv)
1373 c-----------------------------------------------------------------------
1374       include 'epos.inc'
1375       parameter (mxindx=1000,mxre=100,mxma=11,mxmx=6)
1376       common/crema/indx(mxindx),rema(mxre,mxma),rewi(mxre,mxma)
1377      *     ,idmx(mxma,mxmx),icre1(mxre,mxma),icre2(mxre,mxma)
1378       integer ic(2)
1379
1380       ic(1)=000000
1381       ic(2)=000000
1382       if(iabs(id).lt.20)then
1383         if(id.eq.1)then
1384           ic(1)=100000
1385           ic(2)=000000
1386         elseif(id.eq.-1)then
1387           ic(1)=000000
1388           ic(2)=100000
1389         elseif(id.eq.2)then
1390           ic(1)=010000
1391           ic(2)=000000
1392         elseif(id.eq.-2)then
1393           ic(1)=000000
1394           ic(2)=010000
1395         elseif(id.eq.3)then
1396           ic(1)=001000
1397           ic(2)=000000
1398         elseif(id.eq.-3)then
1399           ic(1)=000000
1400           ic(2)=001000
1401         elseif(id.eq.4)then
1402           ic(1)=000100
1403           ic(2)=000000
1404         elseif(id.eq.-4)then
1405           ic(1)=000000
1406           ic(2)=000100
1407         elseif(id.eq.5)then
1408           ic(1)=000010
1409           ic(2)=000000
1410         elseif(id.eq.-5)then
1411           ic(1)=000000
1412           ic(2)=000010
1413         elseif(id.eq.17)then
1414           ic(1)=330000
1415           ic(2)=000000
1416         elseif(id.eq.-17)then
1417           ic(1)=000000
1418           ic(2)=330000
1419         elseif(id.eq.18)then
1420           ic(1)=450000
1421           ic(2)=000000
1422         elseif(id.eq.-18)then
1423           ic(1)=000000
1424           ic(2)=450000
1425         elseif(id.eq.19)then
1426           ic(1)=660000
1427           ic(2)=000000
1428         elseif(id.eq.-19)then
1429           ic(1)=000000
1430           ic(2)=660000
1431         endif
1432         return
1433       endif
1434       if(id.eq.30)then
1435          ic(1)=222000
1436          ic(2)=000000
1437          return
1438       endif
1439       if(iabs(id).lt.1e8)then
1440         ix=iabs(id)/10
1441         if(ix.lt.1.or.ix.gt.mxindx)goto9999
1442         ii=indx(ix)
1443         if(ii.eq.0)goto9998
1444         jj=mod(iabs(id),10)+2
1445         do 27 imx=1,mxmx
1446           do 27 ima=2,mxma
1447             if(iabs(id).eq.idmx(ima,imx))jj=ima
1448  27       continue
1449           if(id.gt.0)then
1450             ic(1)=icre1(ii,jj)
1451             ic(2)=icre2(ii,jj)
1452           else
1453             ic(2)=icre1(ii,jj)
1454             ic(1)=icre2(ii,jj)
1455           endif
1456           if(ic(1).eq.100000.and.ic(2).eq.100000.and.rangen().lt.0.5)
1457      $         then
1458             ic(1)=010000
1459             ic(2)=010000
1460           endif
1461         elseif(mod(id/10**8,10).eq.8)then
1462           ic(1)=mod(id,10**8)/10000*100
1463           ic(2)=mod(id,10**4)*100
1464         else
1465           write(ifch,*)'***** id: ',id
1466           call utstop('idtr4: unrecognized id&')
1467         endif
1468         return
1469
1470  9998   continue
1471         write(ifch,*)'id: ',id
1472         call utstop('idtr4: indx=0.&')
1473
1474  9999   continue
1475         write(ifch,*)'id: ',id
1476         call utstop('idtr4: ix out of range.&')
1477         end
1478
1479 c-----------------------------------------------------------------------
1480       integer function idtra(ic,ier,ires,imix)
1481 c-----------------------------------------------------------------------
1482 c     tranforms from werner-id to paige-id
1483 c         ier .... error message (1) or not (0) in case of problem
1484 c         ires ... dummy variable, not used  any more !!!!
1485 c         imix ... 1 not supported any more
1486 c                  2 010000 010000 -> 110, 001000 000100 -> 110
1487 c                  3 010000 010000 -> 110, 001000 000100 -> 220
1488 c-----------------------------------------------------------------------
1489       include 'epos.inc'
1490       parameter (nidt=51)
1491       integer idt(3,nidt),ic(2)!,icm(2)
1492       data idt/
1493      * 100000,100000, 110   ,100000,010000, 120   ,010000,010000, 220
1494      *,100000,001000, 130   ,010000,001000, 230   ,001000,001000, 330
1495      *,100000,000100, 140   ,010000,000100, 240   ,001000,000100, 340
1496      *,000100,000100, 440
1497      *,100000,000010, 150   ,010000,000010, 250   ,001000,000010, 350
1498      *,000100,000010, 450   ,000010,000010, 550
1499      *,100000,000000,   1   ,010000,000000,   2   ,001000,000000,   3
1500      *,000100,000000,   4   ,000010,000000,   5   ,000001,000000,   6
1501 ccc     *,330000,000000,  17   ,450000,000000,  18   ,660000,000000,  19
1502      *,200000,000000,1100   ,110000,000000,1200   ,020000,000000,2200
1503      *,101000,000000,1300   ,011000,000000,2300   ,002000,000000,3300
1504      *,100100,000000,1400   ,010100,000000,2400   ,001100,000000,3400
1505      *,000200,000000,4400
1506      *,300000,000000,1111   ,210000,000000,1120   ,120000,000000,1220
1507      *,030000,000000,2221   ,201000,000000,1130   ,111000,000000,1230
1508      *,021000,000000,2230   ,102000,000000,1330   ,012000,000000,2330
1509      *,003000,000000,3331   ,200100,000000,1140   ,110100,000000,1240
1510      *,020100,000000,2240   ,101100,000000,1340   ,011100,000000,2340
1511      *,002100,000000,3340   ,100200,000000,1440   ,010200,000000,2440
1512      *,001200,000000,3440   ,000300,000000,4441/
1513
1514       idtra=0
1515       if(ic(1).eq.0.and.ic(2).eq.0)return
1516       i=1
1517       do while(i.le.nidt.and.idtra.eq.0)
1518         if(ic(2).eq.idt(1,i).and.ic(1).eq.idt(2,i))idtra=-idt(3,i)
1519         if(ic(1).eq.idt(1,i).and.ic(2).eq.idt(2,i))idtra=idt(3,i)
1520         i=i+1
1521       enddo
1522       isi=1
1523       if(idtra.ne.0)isi=idtra/iabs(idtra)
1524
1525       jspin=0
1526
1527       if(imix.eq.1)stop'imix=1 no longer supported'
1528       if(imix.eq.2)then
1529       if(idtra.eq.220)idtra=110
1530       if(idtra.eq.330)idtra=110
1531       elseif(imix.eq.3)then
1532       if(idtra.eq.220)idtra=110
1533       if(idtra.eq.330)idtra=220
1534       endif
1535
1536       if(idtra.ne.0)idtra=idtra+jspin*isi
1537
1538       if(idtra.ne.0)return
1539       if(ier.ne.1)return
1540       write(ifch,*)'idtra: ic = ',ic
1541       call utstop('idtra: unknown code&')
1542
1543       entry idtrai(num,id,ier)
1544       idtrai=0
1545       if(iabs(id).eq.20)then
1546       j=5
1547       else
1548       j=0
1549       do 2 i=1,nidt
1550       if(iabs(id).eq.idt(3,i))j=i
1551 2     continue
1552       endif
1553       if(j.ne.0)then
1554       if(id.lt.0)then
1555       idtrai=idt(3-num,j)
1556       else
1557       idtrai=idt(num,j)
1558       endif
1559       return
1560       endif
1561       if(ier.ne.1)return
1562       write(ifch,*)'idtrai: id = ',id
1563       call utstop('idtrai: unknown code&')
1564       end
1565
1566 c-----------------------------------------------------------------------
1567       subroutine idtrb(ib1,ib2,ib3,ib4,jc)
1568 c     id transformation ib -> jc
1569 c-----------------------------------------------------------------------
1570       parameter (nflav=6)
1571       integer jc(nflav,2)
1572       jc(1,1)=ib1/10**4
1573       jc(2,1)=ib2/10**4
1574       jc(3,1)=ib3/10**4
1575       jc(4,1)=ib4/10**4
1576       jc(5,1)=0
1577       jc(6,1)=0
1578       jc(1,2)=mod(ib1,10**4)
1579       jc(2,2)=mod(ib2,10**4)
1580       jc(3,2)=mod(ib3,10**4)
1581       jc(4,2)=mod(ib4,10**4)
1582       jc(5,2)=0
1583       jc(6,2)=0
1584       return
1585       end
1586
1587 c-----------------------------------------------------------------------
1588       subroutine idtrbi(jc,ib1,ib2,ib3,ib4)
1589 c     id transformation jc -> ib
1590 c-----------------------------------------------------------------------
1591       include 'epos.inc'
1592       integer jc(nflav,2)
1593       ib1=jc(1,1)*10**4+jc(1,2)
1594       ib2=jc(2,1)*10**4+jc(2,2)
1595       ib3=jc(3,1)*10**4+jc(3,2)
1596       ib4=jc(4,1)*10**4+jc(4,2)
1597       ib5=jc(5,1)*10**4+jc(5,2)
1598       ib6=jc(6,1)*10**4+jc(6,2)
1599       if(ib5.ne.0.or.ib6.ne.0)then
1600       write(ifch,*)'***** error in idtrbi: bottom or top quarks'
1601       write(ifch,*)'jc:'
1602       write(ifch,*)jc
1603       call utstop('idtrbi: bottom or top quarks&')
1604       endif
1605       return
1606       end
1607
1608 c------------------------------------------------------------------------------
1609       integer function idtrafo(code1,code2,idi)
1610 c------------------------------------------------------------------------------
1611 c.....tranforms id of code1 (=idi) into id of code2 (=idtrafocx)
1612 c.....supported codes:
1613 c.....'nxs' = epos
1614 c.....'pdg' = PDG 1996
1615 c.....'qgs' = QGSJet
1616 c.....'ghe' = Gheisha
1617 c.....'sib' = Sibyll
1618 c.....'cor' = Corsika (GEANT)
1619
1620 C --- ighenex(I)=EPOS CODE CORRESPONDING TO GHEISHA CODE I ---
1621
1622       common /ighnx/ ighenex(35)
1623       data ighenex/
1624      $               10,   11,   -12,    12,   -14,    14,   120,   110,
1625      $             -120,  130,    20,   -20,  -130,  1120, -1120,  1220,
1626      $            -1220, 2130, -2130,  1130,  1230,  2230, -1130, -1230,
1627      $            -2230, 1330,  2330, -1330, -2330,    17,    18,    19,
1628      $            3331, -3331,  30/
1629
1630 C --- DATA STMTS. FOR GEANT/GHEISHA PARTICLE CODE CONVERSIONS ---
1631 C --- KIPART(I)=GHEISHA CODE CORRESPONDING TO GEANT   CODE I ---
1632 C --- IKPART(I)=GEANT   CODE CORRESPONDING TO GHEISHA CODE I ---
1633       DIMENSION        KIPART(48),IKPART(35)
1634       DATA KIPART/
1635      $               1,   3,   4,   2,   5,   6,   8,   7,
1636      $               9,  12,  10,  13,  16,  14,  15,  11,
1637      $              35,  18,  20,  21,  22,  26,  27,  33,
1638      $              17,  19,  23,  24,  25,  28,  29,  34,
1639      $              35,  35,  35,  35,  35,  35,  35,  35,
1640      $              35,  35,  35,  35,  30,  31,  32,  35/
1641
1642       DATA IKPART/
1643      $               1,   4,   2,   3,   5,   6,   8,   7,
1644      $               9,  11,  16,  10,  12,  14,  15,  13,
1645      $              25,  18,  26,  19,  20,  21,  27,  28,
1646      $              29,  22,  23,  30,  31,  45,  46,  47,
1647      $              24,  32,  48/
1648 c-------------------------------------------------------------------------------
1649
1650       character*3 code1,code2
1651       parameter (ncode=5,nidt=334)
1652       integer idt(ncode,nidt)
1653       double precision drangen,dummy
1654
1655       data ((idt(i,j),i=1,ncode),j= 1,18)/
1656      *          1,2,99,99,99             !u quark
1657      *     ,    2,1,99,99,99             !d
1658      *     ,    3,3,99,99,99             !s
1659      *     ,    4,4,99,99,99             !c
1660      *     ,    5,5,99,99,99             !b
1661      *     ,    6,6,99,99,99             !t
1662      *     ,   10,22,99,1,1              !gamma
1663      *     ,   9 ,21,99,99,99            !gluon
1664      *     ,   12,11,11,4,3              !e-
1665      *     ,  -12,-11,-11,3,2            !e+
1666      *     ,   11,12,99,2,15             !nu_e-
1667      *     ,  -11,-12,99,-2,16           !nu_e+
1668      *     ,   14,13,99,6,5              !mu-
1669      *     ,  -14,-13,99,5,4             !mu+
1670      *     ,   13,14,99,2,17             !nu_mu-
1671      *     ,  -13,-14,99,-2,18           !nu_mu+
1672      *     ,   16,15,99,99,19            !tau-
1673      *     ,   15,16,99,99,20 /          !nu_tau-
1674       data ((idt(i,j),i=1,ncode),j= 19,40)/
1675      *        110,111,0,8,6              !pi0
1676      *     ,  120,211,1,7,7              !pi+
1677      *     , -120,-211,-1,9,8            !pi-
1678      *     ,  220,221,10,99,23           !eta
1679      *     ,  130,321,4,10,9             !k+
1680      *     , -130,-321,-4,13,10          !k-
1681      *     ,  230,311,5,11,21            !k0
1682      *     , -230,-311,-5,12,22          !k0b
1683      *     ,   20,310,5,11,12            !kshort
1684      *     ,  -20,-310,-5,12,11          !klong
1685      *     ,  330,331,99,99,24           !etaprime
1686      *     ,  111,113,99,99,27           !rho0
1687      *     ,  121,213,99,99,25           !rho+
1688      *     , -121,-213,99,99,26          !rho-
1689      *     ,  221,223,99,99,32           !omega
1690      *     ,  131,323,99,99,28           !k*+
1691      *     , -131,-323,99,99,29          !k*-
1692      *     ,  231,313,99,99,30           !k*0
1693      *     , -231,-313,99,99,31          !k*0b
1694      *     ,  331,333,99,99,33           !phi
1695      $     , -140,421,8,99,99            !D0(1.864)
1696      $     ,  240,-411,7,99,99 /         !D(1.869)-
1697       data ((idt(i,j),i=1,ncode),j= 41,59)/
1698      *       1120,2212,2,14,13           !proton
1699      *     , 1220,2112,3,16,14           !neutron
1700      *     , 2130,3122,6,18,39           !lambda
1701      *     , 1130,3222,99,20,34          !sigma+
1702      *     , 1230,3212,99,21,35          !sigma0
1703      *     , 2230,3112,99,22,36          !sigma-
1704      *     , 1330,3322,99,26,37          !xi0
1705      *     , 2330,3312,99,27,38          !xi-
1706      *     , 1111,2224,99,99,40          !delta++
1707      *     , 1121,2214,99,99,41          !delta+
1708      *     , 1221,2114,99,99,42          !delta0
1709      *     , 2221,1114,99,99,43          !delta-
1710      *     , 1131,3224,99,99,44          !sigma*+
1711      *     , 1231,3214,99,99,45          !sigma*0
1712      *     , 2231,3114,99,99,46          !sigma*-
1713      *     , 1331, 3324,99,99,47         !xi*0
1714      *     , 2331, 3314,99,99,48         !xi*-
1715      *     , 3331, 3334,99,33,49         !omega-
1716      $     , 2140, 4122,9,99,99   /      !LambdaC(2.285)+
1717       data ((idt(i,j),i=1,ncode),j= 60,64)/
1718      $      17,99,99,30,1002             !  Deuteron
1719      $     ,18,99,99,31,1003             !  Triton
1720      $     ,19,99,99,32,1004             !  Alpha
1721      $     ,0,99,99,0,0                  !  Air
1722      *     ,99,99,99,99,99 /             !  unknown
1723       data ((idt(i,j),i=1,ncode),j= 65,79)/
1724      $       1112,32224,99,99,99         !  Delta(1600)++
1725      $     , 1112, 2222,99,99,99         !  Delta(1620)++
1726      $     , 1113,12224,99,99,99         !  Delta(1700)++
1727      $     , 1114,12222,99,99,99         !  Delta(1900)++
1728      $     , 1114, 2226,99,99,99         !  Delta(1905)++
1729      $     , 1114,22222,99,99,99         !  Delta(1910)++
1730      $     , 1114,22224,99,99,99         !  Delta(1920)++
1731      $     , 1114,12226,99,99,99         !  Delta(1930)++
1732      $     , 1114, 2228,99,99,99         !  Delta(1950)++
1733      $     , 2222,31114,99,99,99         !  Delta(1600)-
1734      $     , 2222, 1112,99,99,99         !  Delta(1620)-
1735      $     , 2223,11114,99,99,99         !  Delta(1700)-
1736      $     , 2224,11112,99,99,99         !  Delta(1900)-
1737      $     , 2224, 1116,99,99,99         !  Delta(1905)-
1738      $     , 2224,21112,99,99,99   /     !  Delta(1910)-
1739       data ((idt(i,j),i=1,ncode),j= 80,94)/
1740      $      2224,21114,99,99,99          !  Delta(1920)-
1741      $     ,2224,11116,99,99,99          !  Delta(1930)-
1742      $     ,2224, 1118,99,99,99          !  Delta(1950)-
1743      $     ,1122,12212,99,99,99          !  N(1440)+
1744      $     ,1123, 2124,99,99,99          !  N(1520)+
1745      $     ,1123,22212,99,99,99          !  N(1535)+
1746      $     ,1124,32214,99,99,99          !  Delta(1600)+
1747      $     ,1124, 2122,99,99,99          !  Delta(1620)+
1748      $     ,1125,32212,99,99,99          !  N(1650)+
1749      $     ,1125, 2216,99,99,99          !  N(1675)+
1750      $     ,1125,12216,99,99,99          !  N(1680)+
1751      $     ,1126,12214,99,99,99          !  Delta(1700)+
1752      $     ,1127,22124,99,99,99          !  N(1700)+
1753      $     ,1127,42212,99,99,99          !  N(1710)+
1754      $     ,1127,32124,99,99,99   /      !  N(1720)+
1755       data ((idt(i,j),i=1,ncode),j= 95,109)/
1756      $      1128,12122,99,99,99          !  Delta(1900)+
1757      $     ,1128, 2126,99,99,99          !  Delta(1905)+
1758      $     ,1128,22122,99,99,99          !  Delta(1910)+
1759      $     ,1128,22214,99,99,99          !  Delta(1920)+
1760      $     ,1128,12126,99,99,99          !  Delta(1930)+
1761      $     ,1128, 2218,99,99,99          !  Delta(1950)+
1762      $     ,1222,12112,99,99,99          !  N(1440)0
1763      $     ,1223, 1214,99,99,99          !  N(1520)0
1764      $     ,1223,22112,99,99,99          !  N(1535)0
1765      $     ,1224,32114,99,99,99          !  Delta(1600)0
1766      $     ,1224, 1212,99,99,99          !  Delta(1620)0
1767      $     ,1225,32112,99,99,99          !  N(1650)0
1768      $     ,1225, 2116,99,99,99          !  N(1675)0
1769      $     ,1225,12116,99,99,99          !  N(1680)0
1770      $     ,1226,12114,99,99,99   /      !  Delta(1700)0
1771       data ((idt(i,j),i=1,ncode),j= 110,124)/
1772      $      1227,21214,99,99,99          !  N(1700)0
1773      $     ,1227,42112,99,99,99          !  N(1710)0
1774      $     ,1227,31214,99,99,99          !  N(1720)0
1775      $     ,1228,11212,99,99,99          !  Delta(1900)0
1776      $     ,1228, 1216,99,99,99          !  Delta(1905)0
1777      $     ,1228,21212,99,99,99          !  Delta(1910)0
1778      $     ,1228,22114,99,99,99          !  Delta(1920)0
1779      $     ,1228,11216,99,99,99          !  Delta(1930)0
1780      $     ,1228, 2118,99,99,99          !  Delta(1950)0
1781      $     ,1233,13122,99,99,99          !  Lambda(1405)0
1782      $     ,1234, 3124,99,99,99          !  Lambda(1520)0
1783      $     ,1235,23122,99,99,99          !  Lambda(1600)0
1784      $     ,1235,33122,99,99,99          !  Lambda(1670)0
1785      $     ,1235,13124,99,99,99          !  Lambda(1690)0
1786      $     ,1236,13212,99,99,99  /       !  Sigma(1660)0
1787       data ((idt(i,j),i=1,ncode),j= 125,139)/
1788      $      1236,13214,99,99,99          !  Sigma(1670)0
1789      $     ,1237,23212,99,99,99          !  Sigma(1750)0
1790      $     ,1237, 3216,99,99,99          !  Sigma(1775)0
1791      $     ,1238,43122,99,99,99          !  Lambda(1800)0
1792      $     ,1238,53122,99,99,99          !  Lambda(1810)0
1793      $     ,1238, 3126,99,99,99          !  Lambda(1820)0
1794      $     ,1238,13126,99,99,99          !  Lambda(1830)0
1795      $     ,1238,23124,99,99,99          !  Lambda(1890)0
1796      $     ,1239,13216,99,99,99          !  Sigma(1915)0
1797      $     ,1239,23214,99,99,99          !  Sigma(1940)0
1798      $     ,1132,13222,99,99,99          !  Sigma(1660)+
1799      $     ,1132,13224,99,99,99          !  Sigma(1670)+
1800      $     ,1133,23222,99,99,99          !  Sigma(1750)+
1801      $     ,1133,3226,99,99,99           !  Sigma(1775)+
1802      $     ,1134,13226,99,99,99   /      !  Sigma(1915)+
1803       data ((idt(i,j),i=1,ncode),j= 140,146)/
1804      $      1134,23224,99,99,99          !  Sigma(1940)+
1805      $     ,2232,13112,99,99,99          !  Sigma(1660)-
1806      $     ,2232,13114,99,99,99          !  Sigma(1670)-
1807      $     ,2233,23112,99,99,99          !  Sigma(1750)-
1808      $     ,2233,3116,99,99,99           !  Sigma(1775)-
1809      $     ,2234,13116,99,99,99          !  Sigma(1915)-
1810      $     ,2234,23114,99,99,99  /       !  Sigma(1940)-
1811       data ((idt(i,j),i=1,ncode),j= 147,159)/
1812      $      5,7,99,99,99                 !  quark b'
1813      $     ,6,8,99,99,99                 !  quark t'
1814      $     ,16,17,99,99,99               !  lepton tau'
1815      $     ,15,18,99,99,99               !  lepton nu' tau
1816      $     ,90,23,99,99,99               !  Z0
1817      $     ,80,24,99,99,99               !  W+
1818      $     ,81,25,99,99,99               !  h0
1819      $     ,85,32,99,99,99               !  Z'0
1820      $     ,86,33,99,99,99               !  Z''0
1821      $     ,87,34,99,99,99               !  W'+
1822      $     ,82,35,99,99,99               !  H0
1823      $     ,83,36,99,99,99               !  A0
1824      $     ,84,37,99,99,99   /           !  H+
1825       data ((idt(i,j),i=1,ncode),j= 160,184)/
1826      $      1200,2101,99,99,99           !  diquark ud_0
1827      $     ,2300,3101,99,99,99           !  diquark sd_0
1828      $     ,1300,3201,99,99,99           !  diquark su_0
1829      $     ,2400,4101,99,99,99           !  diquark cd_0
1830      $     ,1400,4201,99,99,99           !  diquark cu_0
1831      $     ,3400,4301,99,99,99           !  diquark cs_0
1832      $     ,2500,5101,99,99,99           !  diquark bd_0
1833      $     ,1500,5201,99,99,99           !  diquark bu_0
1834      $     ,3500,5301,99,99,99           !  diquark bs_0
1835      $     ,4500,5401,99,99,99           !  diquark bc_0
1836      $     ,2200,1103,99,99,99           !  diquark dd_1
1837      $     ,1200,2103,99,99,99           !  diquark ud_1
1838      $     ,1100,2203,99,99,99           !  diquark uu_1
1839      $     ,2300,3103,99,99,99           !  diquark sd_1
1840      $     ,1300,3203,99,99,99           !  diquark su_1
1841      $     ,3300,3303,99,99,99           !  diquark ss_1
1842      $     ,2400,4103,99,99,99           !  diquark cd_1
1843      $     ,1400,4203,99,99,99           !  diquark cu_1
1844      $     ,3400,4303,99,99,99           !  diquark cs_1
1845      $     ,4400,4403,99,99,99           !  diquark cc_1
1846      $     ,2500,5103,99,99,99           !  diquark bd_1
1847      $     ,1500,5203,99,99,99           !  diquark bu_1
1848      $     ,3500,5303,99,99,99           !  diquark bs_1
1849      $     ,4500,5403,99,99,99           !  diquark bc_1
1850      $     ,5500,5503,99,99,99 /         !  diquark bb_1
1851       data ((idt(i,j),i=1,ncode),j= 185,188)/
1852      $      800000091,91,99,99,99      !  parton system in cluster fragmentation  (pythia)
1853      $     ,800000092,92,99,99,99      !  parton system in string fragmentation  (pythia)
1854      $     ,800000093,93,99,99,99      !  parton system in independent system  (pythia)
1855      $     ,800000094,94,99,99,99 /    !  CMshower (pythia)
1856       data ((idt(i,j),i=1,ncode),j= 189,208)/
1857      $      -340,431,99,99,99            !  Ds+
1858      $     ,340,-431,99,99,99            !  Ds-
1859      $     ,-241,413,99,99,99            !  D*+
1860      $     ,241,-413,99,99,99            !  D*-
1861      $     ,-141,423,99,99,99            !  D*0
1862      $     ,141,-423,99,99,99            !  D*0b
1863      $     ,-341,433,99,99,99            !  Ds*+
1864      $     ,341,-433,99,99,99            !  Ds*-
1865      $     ,250,511,99,99,99             !  B0
1866      $     ,150,521,99,99,99             !  B+
1867      $     ,350,531,99,99,99             !  B0s+
1868      $     ,450,541,99,99,99             !  Bc+
1869      $     ,251,513,99,99,99             !  B*0
1870      $     ,151,523,99,99,99             !  B*+
1871      $     ,351,533,99,99,99             !  B*0s+
1872      $     ,451,543,99,99,99             !  B*c+
1873      $     ,440,441,99,99,99             !  etac
1874      $     ,441,443,99,99,99             !  J/psi
1875      $     ,550,551,99,99,99             !  etab
1876      $     ,551,553,99,99,99   /         !  Upsilon
1877       data ((idt(i,j),i=1,ncode),j= 209,264)/
1878      $      2240,4112,99,99,99           !  sigmac0
1879      $     ,1240,4212,99,99,99           !  sigmac+
1880      $     ,1140,4222,99,99,99           !  sigmac++
1881      $     ,2241,4114,99,99,99           !  sigma*c0
1882      $     ,1241,4214,99,99,99           !  sigma*c+
1883      $     ,1141,4224,99,99,99           !  sigma*c++
1884      $     ,3240,4132,99,99,99           !  Xic0
1885      $     ,2340,4312,99,99,99           !  Xi'c0
1886      $     ,3140,4232,99,99,99           !  Xic+
1887      $     ,1340,4322,99,99,99           !  Xi'c+
1888      $     ,3340,4332,99,99,99           !  omegac0
1889      $     ,2341,4314,99,99,99           !  Xi*c0
1890      $     ,1341,4324,99,99,99           !  Xi*c+
1891      $     ,3341,4334,99,99,99           !  omega*c0
1892      $     ,2440,4412,99,99,99           !  dcc
1893      $     ,2441,4414,99,99,99           !  dcc*
1894      $     ,1440,4422,99,99,99           !  ucc
1895      $     ,1441,4424,99,99,99           !  ucc*
1896      $     ,3440,4432,99,99,99           !  scc
1897      $     ,3441,4434,99,99,99           !  scc*
1898      $     ,4441,4444,99,99,99           !  ccc*
1899      $     ,2250,5112,99,99,99           !  sigmab-
1900      $     ,2150,5122,99,99,99           !  lambdab0
1901      $     ,3250,5132,99,99,99           !  sdb
1902      $     ,4250,5142,99,99,99           !  cdb
1903      $     ,1250,5212,99,99,99           !  sigmab0
1904      $     ,1150,5222,99,99,99           !  sigmab+
1905      $     ,3150,5232,99,99,99           !  sub
1906      $     ,4150,5242,99,99,99           !  cub
1907      $     ,2350,5312,99,99,99           !  dsb
1908      $     ,1350,5322,99,99,99           !  usb
1909      $     ,3350,5332,99,99,99           !  ssb
1910      $     ,4350,5342,99,99,99           !  csb
1911      $     ,2450,5412,99,99,99           !  dcb
1912      $     ,1450,5422,99,99,99           !  ucb
1913      $     ,3450,5432,99,99,99           !  scb
1914      $     ,4450,5442,99,99,99           !  ccb
1915      $     ,2550,5512,99,99,99           !  dbb
1916      $     ,1550,5522,99,99,99           !  ubb
1917      $     ,3550,5532,99,99,99           !  sbb
1918      $     ,3550,5542,99,99,99           !  scb
1919      $     ,2251,5114,99,99,99           !  sigma*b-
1920      $     ,1251,5214,99,99,99           !  sigma*b0
1921      $     ,1151,5224,99,99,99           !  sigma*b+
1922      $     ,2351,5314,99,99,99           !  dsb*
1923      $     ,1351,5324,99,99,99           !  usb*
1924      $     ,3351,5334,99,99,99           !  ssb*
1925      $     ,2451,5414,99,99,99           !  dcb*
1926      $     ,1451,5424,99,99,99           !  ucb*
1927      $     ,3451,5434,99,99,99           !  scb*
1928      $     ,4451,5444,99,99,99           !  ccb*
1929      $     ,2551,5514,99,99,99           !  dbb*
1930      $     ,1551,5524,99,99,99           !  ubb*
1931      $     ,3551,5534,99,99,99           !  sbb*
1932      $     ,4551,5544,99,99,99           !  cbb*
1933      $     ,5551,5554,99,99,99  /        !  bbb*
1934       data ((idt(i,j),i=1,ncode),j= 265,295)/
1935      $      123,10213,99,99,99           !  b1
1936      $     ,122,10211,99,99,99           !  a0+
1937      $     ,233,10313,99,99,99           !  K0_1
1938      $     ,232,10311,99,99,99           !  K*0_1
1939      $     ,133,10323,99,99,99           !  K+_1
1940      $     ,132,10321,99,99,99           !  K*+_1
1941      $     ,143,10423,99,99,99           !  D0_1
1942      $     ,132,10421,99,99,99           !  D*0_1
1943      $     ,243,10413,99,99,99           !  D+_1
1944      $     ,242,10411,99,99,99           !  D*+_1
1945      $     ,343,10433,99,99,99           !  D+s_1
1946      $     ,342,10431,99,99,99           !  D*0s+_1
1947      $     ,223,10113,99,99,99           !  b_10
1948      $     ,222,10111,99,99,99           !  a_00
1949      $     ,113,10223,99,99,99           !  h_10
1950      $     ,112,10221,99,99,99           !  f_00
1951      $     ,333,10333,99,99,99           !  h'_10
1952      $     ,332,10331,99,99,99           !  f'_00
1953      $     ,443,10443,99,99,99           !  h_1c0
1954      $     ,442,10441,99,99,99           !  Xi_0c0
1955      $     ,444,10443,99,99,99           !  psi'
1956      $     ,253,10513,99,99,99           !  db_10
1957      $     ,252,10511,99,99,99           !  db*_00
1958      $     ,153,10523,99,99,99           !  ub_10
1959      $     ,152,10521,99,99,99           !  ub*_00
1960      $     ,353,10533,99,99,99           !  sb_10
1961      $     ,352,10531,99,99,99           !  sb*_00
1962      $     ,453,10543,99,99,99           !  cb_10
1963      $     ,452,10541,99,99,99           !  cb*_00
1964      $     ,553,10553,99,99,99           !  Upsilon'
1965      $     ,552,10551,99,99,99  /        !  Upsilon'*
1966       data ((idt(i,j),i=1,ncode),j= 296,325)/
1967      $      124,20213,99,99,99           !  a_1+
1968      $     ,125,215,99,99,99             !  a_2+
1969      $     ,234,20313,99,99,99           !  K*0_1
1970      $     ,235,315,99,99,99             !  K*0_2
1971      $     ,134,20323,99,99,99           !  K*+_1
1972      $     ,135,325,99,99,99             !  K*+_2
1973      $     ,144,20423,99,99,99           !  D*_10
1974      $     ,135,425,99,99,99             !  D*_20
1975      $     ,244,20413,99,99,99           !  D*_1+
1976      $     ,245,415,99,99,99             !  D*_2+
1977      $     ,344,20433,99,99,99           !  D*_1s+
1978      $     ,345,435,99,99,99             !  D*_2s+
1979      $     ,224,20113,99,99,99           !  a_10
1980      $     ,225,115,99,99,99             !  a_20
1981      $     ,114,20223,99,99,99           !  f_10
1982      $     ,115,225,99,99,99             !  f_20
1983      $     ,334,20333,99,99,99           !  f'_10
1984      $     ,335,335,99,99,99             !  f'_20
1985      $     ,444,20443,99,99,99           !  Xi_1c0
1986      $     ,445,445,99,99,99             !  Xi_2c0
1987      $     ,254,20513,99,99,99           !  db*_10
1988      $     ,255,515,99,99,99             !  db*_20
1989      $     ,154,20523,99,99,99           !  ub*_10
1990      $     ,155,525,99,99,99             !  ub*_20
1991      $     ,354,20533,99,99,99           !  sb*_10
1992      $     ,355,535,99,99,99             !  sb*_20
1993      $     ,454,20543,99,99,99           !  cb*_10
1994      $     ,455,545,99,99,99             !  cb*_20
1995      $     ,554,20553,99,99,99           !  bb*_10
1996      $     ,555,555,99,99,99    /        !  bb*_20
1997       data ((idt(i,j),i=1,ncode),j= 326,nidt)/
1998      $      11099,9900110,99,99,99       !  diff pi0 state
1999      $     ,12099,9900210,99,99,99       !  diff pi+ state
2000      $     ,22099,9900220,99,99,99       !  diff omega state
2001      $     ,33099,9900330,99,99,99       !  diff phi state
2002      $     ,44099,9900440,99,99,99       !  diff J/psi state
2003      $     ,112099,9902210,99,99,99      !  diff proton state
2004      $     ,122099,9902110,99,99,99      !  diff neutron state
2005      $     ,800000110,110,99,99,99       !  Reggeon
2006      $     ,800000990,990,99,99,99 /      !  Pomeron
2007
2008
2009
2010
2011       nidtmx=64
2012       id1=idi
2013       if(code1.eq.'nxs')then
2014         i=1
2015       elseif(code1.eq.'pdg')then
2016         i=2
2017       elseif(code1.eq.'qgs')then
2018         i=3
2019       elseif(code1.eq.'ghe')then
2020         id1=ighenex(id1)
2021         i=1
2022       elseif(code1.eq.'sib')then
2023         i=5
2024       elseif(code1.eq.'cor')then
2025         id1=kipart(id1)
2026         id1=ighenex(id1)
2027         i=1
2028       else
2029         stop "unknown code in idtrafo"
2030       endif
2031       if(code2.eq.'nxs')then
2032         j=1
2033         ji=j
2034         if(i.eq.5.and.id1.gt.1004)then               !nucleus from Sibyll
2035           idtrafo=(id1-1000)*100
2036           return
2037         elseif(id1.eq.130.and.i.eq.2)then
2038           idtrafo=-20
2039           return
2040         endif
2041         if(i.eq.2) nidtmx=nidt
2042       elseif(code2.eq.'pdg')then
2043         j=2
2044         ji=j
2045         if(id1.eq.-20.and.i.eq.1)then
2046           idtrafo=130
2047           return
2048         endif
2049         if(i.eq.1) nidtmx=nidt
2050       elseif(code2.eq.'qgs')then
2051         j=3
2052         ji=j
2053       elseif(code2.eq.'ghe')then
2054         j=4
2055         ji=j
2056       elseif(code2.eq.'sib')then
2057         j=5
2058         ji=j
2059       elseif(code2.eq.'cor')then
2060         j=4
2061         ji=6
2062       else
2063         stop "unknown code in idtrafo"
2064       endif
2065       iad1=abs(id1)
2066       isi=sign(1,id1)
2067
2068       do n=1,nidtmx
2069         if(iad1.eq.abs(idt(i,n)))then
2070           m=1
2071           do while(abs(idt(i,n+m)).eq.iad1)
2072             m=m+1
2073           enddo
2074           mm=0
2075           if(m.gt.1)then
2076             if(m.eq.2.and.idt(i,n)*idt(i,n+1).lt.0)then
2077               if(id1.eq.idt(i,n+1))mm=1
2078               isi=1
2079             else
2080               mm=int(drangen(dummy)*dble(m))
2081             endif
2082           endif
2083           idtrafo=idt(j,n+mm)*isi
2084           if(abs(idtrafo).eq.99)call utstop('New particle not allowed ')
2085           if(idtrafo.lt.0.and.j.eq.4)then                      !gheisha  id always >0
2086             iadtr=abs(idtrafo)
2087             if(iadtr.ge.20.and.iadtr.le.22)then
2088               idtrafo=iadtr+3
2089             elseif(iadtr.eq.26.or.iadtr.eq.27)then
2090               idtrafo=iadtr+2
2091             elseif(iadtr.ge.14)then
2092               idtrafo=iadtr+1
2093             else
2094               idtrafo=iadtr
2095             endif
2096           endif
2097           if(ji.eq.6)idtrafo=ikpart(idtrafo)
2098           return
2099         end if
2100       enddo
2101
2102       print *, 'particle:',code1,'->', code2,id1
2103       stop'idtrafo: nothing found'
2104
2105       end
2106