string.f

c $Id: string.f,v 1.18 2003/05/02 11:06:56 weber Exp $
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
      subroutine stringdec(ityp,iz2,smass,part,ident2,npart)
c
cinput smass   : Stringmass
cinput ityp    : Particle ID
cinput iz2     : Isospin$_3\cdot 2$
c
coutput part   : 4-momenta, 4-position, masses (array)
coutput ident2 : ityp, iz2 (array)
coutput npart  : number of outgoing particles
c
c     This subroutine performs string fragmentation.
c
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      implicit real*8 (a-h,o-z)
      implicit integer (i-n)

      PARAMETER(MXPTCL=200)
      COMMON/PARTCL/ PPTCL(9,MXPTCL),nptcl,IDENT(MXPTCL),IDCAY(MXPTCL)

      include 'comstr.f'

      dimension part(9,mxptcl)
      dimension ident2(2,mxptcl)


c we call the translation routine
        call ityp2id(ityp,iz2,ifa,ifb)

        goto 1
cspl... string fragmentation called with quark id's and energy as arguments
         entry qstring(ifanew,ifbnew,smass,part,ident2,npart)
         ifa=ifanew
         ifb=ifbnew
 1      continue


        smem=smass
c here we call the fragmentation routine. the produced hadrons and their
c properties are returned via the pptcl- and ident-array in the common-block
       call string(ifa,ifb,smass)
c here the array pptcl has been filled with nptcl entries (1->nptcl)
c now we translate to uqmd and shift the pptcl- and ident-info to the
c corresponding part- and ident2-arrays of the newpart-common-block:
       npart=nptcl
       do 2 i=1,nptcl
        call id2ityp(ident(i),pptcl(5,i),itypout,iz2out)
        ident2(1,i)=itypout
        ident2(2,i)=iz2out
        smem=smem-pptcl(4,i)
        do 3 j=1,9
         part(j,i)=pptcl(j,i)
 3       continue
 2      continue
c check for energy conservation:
ctp060926       if(abs(smem).gt.1.0D-5)then
ctp060926        write(*,*)'! stringdec: energy difference=',smem
ctp060926        write(*,*)'ifa,ifb,smass=',ifa,ifb,smass
ctp060926       endif
       return

       end

cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
      subroutine strini
c
c output     : via common blocks
c
c     {\tt strini} calculates mixing angles for the meson-multipletts
c
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      implicit real*8 (a-h,o-z)
      implicit integer (i-n)

      include 'options.f'
        include 'comres.f'
      include 'comstr.f'
        real*8 m3
        real*8 massit
        integer jit
c mixing angles of meson multiplets according to flavor SU(3) quark model:
cspl-0795 these parameters assign the pure u/ubar,d/dbar,s/sbar states
c  (e.g. 110,220,330) to the physical particles according to the su(3)
c  quark model. The flavor mixing angles are chosen according to quadratic
c  Gell-Mann-Okubo mass formula (Review of Particle Properties,
c  Phys. Rev D50 (1994) 1319). For the scalar mesons this formula is not
c  applicable. We assume an ideal mixing angle (tan(theta)=1/sqrt(2)).
c
c pseudoscalar: theta=-10 deg
c vector      : theta= 39 deg
c pseudovector: theta= 51 deg
c tensor      : theta= 28 deg
c
      real*8 mixang(njspin)
c ideal mixing angles assumed for the last three multiplets
      data mixang/-10d0,39d0,35.3d0,51d0,28d0,35.3d0,35.3d0,35.3d0/
c

      pi = 4d0*datan2(1d0,1d0)

c calculate 'singlet shift probabilities', e.g. a 11 (u-ubar) state
c can be changed to a 22 (d-dbar) or 33 (s-sbar) state with a certain
c probability. THEN they can be identified with physical hadrons!
      do 3 i=1,njspin
         mixang(i)=mixang(i)/36d1*2d0*3.1416d0
         PMIX1S(1,i)=(dcos(mixang(i))/sqrt(6d0)
     &        -dsin(mixang(i))/sqrt(3d0))**2
         PMIX1S(2,i)=PMIX1S(1,i)
         PMIX1S(3,i)=(-(dcos(mixang(i))*2d0/dsqrt(6d0))
     &        -dsin(mixang(i))/dsqrt(3d0))**2
         PMIX2S(1,i)=0.5d0
         PMIX2S(2,i)=0.5d0
         PMIX2S(3,i)=1d0

ce calculate probabilities of the meson multipletts
ce according to parm=(spin degeneracy)/(average mass) *ctp(50 ff.)
           parm(i)=0d0
           m3=0d0
           do 102 j=0,3
             itp=mlt2it(4*(i-1)+j+1)
             m3=m3+massit(itp)
             jpc=jit(itp)/2
102           continue
           parm(i)=parm(i)+(2*jpc+1)/m3*4*CTParam(49+i)

c the mixing-angles are the same for 'string' and 'cluster':
         do 2 k=1,3
            PMIX1C(k,i)=PMIX1S(k,i)
            PMIX2C(k,i)=PMIX2S(k,i)
c            write(6,*)'#',pmix1c(k,i)
 2       continue
 3    continue

      return
      end

cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
      SUBROUTINE GAUSPT(PT0,SIGQT)
c
cinput   sigqt  : Width of Gaussian
c
coutput  pt0    : transverse momentum
c
C     generate pt with Gaussian
c     distribution $\propto pt \exp(-pt^2/sigqt^2)$
c
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      implicit none
      real*8 pt0,sigqt,rnd,ranf

      RND=ranf(0)
      PT0=SIGQT*SQRT(-DLOG(1.d0-RND))
      RETURN
      END


cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
      SUBROUTINE FLAVOR(ID,IFL1,IFL2,IFL3,JSPIN)
c
cinput ID     : quarkcode
c
coutput ifl1  : single quarks id
coutput ifl2  : single quarks id
coutput ifl3  : single quarks id
coutput jspin : spin id
c
C          THIS SUBROUTINE UNPACKS THE IDENT CODE ID=+/-IJKL
c
C          -MESONS:
C          I=0, J<=K, +/- IS SIGN FOR J,
C          ID=110 FOR PI0, ID=220 FOR ETA, ETC.
c
C          -BARYONS:
C          I<=J<=K IN GENERAL,
C          J<I<K FOR SECOND STATE ANTISYMMETRIC IN (I,J), EG. L = 2130
c
C          -DIQUARKS:
C          ID=+/-IJ00, I<J FOR DIQUARK COMPOSED OF I,J.
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      implicit real*8 (a-h,o-z)
      implicit integer (i-n)

      IDABS=IABS(ID)
c extract the single (anti-)quark ids from the hadron id:
      I=IDABS/1000
      J=MOD(IDABS/100,10)
      K=MOD(IDABS/10,10)
      JSPIN=MOD(IDABS,10)
c diquarks:
      IF(ID.NE.0.AND.MOD(ID,100).EQ.0) GO TO 300
c quarks oder so:
      IF(J.EQ.0) GO TO 200
c mesonen:
      IF(I.EQ.0) GO TO 100

C..BARYONS:
C..ONLY X,Y BARYONS ARE QQX, QQY, Q=U,D,S.
      IFL1=ISIGN(I,ID)
      IFL2=ISIGN(J,ID)