lueevt.f

Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file lueevt.f

    
C*********************************************************************  
    
      SUBROUTINE lueevt(KFL,ECM)    
    
C...Purpose: to handle the generation of an e+e- annihilation jet event.    
      IMPLICIT DOUBLE PRECISION(d)  
      common/lujets/n,k(9000,5),p(9000,5),v(9000,5)
      SAVE /lujets/ 
      common/ludat1/mstu(200),paru(200),mstj(200),parj(200) 
      SAVE /ludat1/ 
      common/ludat2/kchg(500,3),pmas(500,4),parf(2000),vckm(4,4)    
      SAVE /ludat2/ 
    
C...Check input parameters. 
      IF(mstu(12).GE.1) CALL lulist(0)  
      IF(kfl.LT.0.OR.kfl.GT.8) THEN 
        CALL luerrm(16,'(LUEEVT:) called with unknown flavour code')    
        IF(mstu(21).GE.1) RETURN    
      ENDIF 
      IF(kfl.LE.5) ecmmin=parj(127)+2.02*parf(100+max(1,kfl))   
      IF(kfl.GE.6) ecmmin=parj(127)+2.02*pmas(kfl,1)    
      IF(ecm.LT.ecmmin) THEN    
        CALL luerrm(16,'(LUEEVT:) called with too small CM energy') 
        IF(mstu(21).GE.1) RETURN    
      ENDIF 
    
C...Check consistency of MSTJ options set.  
      IF(mstj(109).EQ.2.AND.mstj(110).NE.1) THEN    
        CALL luerrm(6,  
     &  '(LUEEVT:) MSTJ(109) value requires MSTJ(110) = 1') 
        mstj(110)=1 
      ENDIF 
      IF(mstj(109).EQ.2.AND.mstj(111).NE.0) THEN    
        CALL luerrm(6,  
     &  '(LUEEVT:) MSTJ(109) value requires MSTJ(111) = 0') 
        mstj(111)=0 
      ENDIF 
    
C...Initialize alpha_strong and total cross-section.    
      mstu(111)=mstj(108)   
      IF(mstj(108).EQ.2.AND.(mstj(101).EQ.0.OR.mstj(101).EQ.1)) 
     &mstu(111)=1   
      paru(112)=parj(121)   
      IF(mstu(111).EQ.2) paru(112)=parj(122)    
      IF(mstj(116).GT.0.AND.(mstj(116).GE.2.OR.abs(ecm-parj(151)).GE.   
     &parj(139).OR.10*mstj(102)+kfl.NE.mstj(119))) CALL luxtot(kfl,ecm, 
     &xtot) 
      IF(mstj(116).GE.3) mstj(116)=1    
    
C...Add initial e+e- to event record (documentation only).  
      ntry=0    
  100 ntry=ntry+1   
      IF(ntry.GT.100) THEN  
        CALL luerrm(14,'(LUEEVT:) caught in an infinite loop')  
        RETURN  
      ENDIF 
      nc=0  
      IF(mstj(115).GE.2) THEN   
        nc=nc+2 
        CALL lu1ent(nc-1,11,0.5*ecm,0.,0.)  
        k(nc-1,1)=21    
        CALL lu1ent(nc,-11,0.5*ecm,paru(1),0.)  
        k(nc,1)=21  
      ENDIF 
    
C...Radiative photon (in initial state).    
      mk=0  
      ecmc=ecm  
      IF(mstj(107).GE.1.AND.mstj(116).GE.1) CALL luradk(ecm,mk,pak, 
     &thek,phik,alpk)   
      IF(mk.EQ.1) ecmc=sqrt(ecm*(ecm-2.*pak))   
      IF(mstj(115).GE.1.AND.mk.EQ.1) THEN   
        nc=nc+1 
        CALL lu1ent(nc,22,pak,thek,phik)    
        k(nc,3)=min(mstj(115)/2,1)  
      ENDIF 
    
C...Virtual exchange boson (gamma or Z0).   
      IF(mstj(115).GE.3) THEN   
        nc=nc+1 
        kf=22   
        IF(mstj(102).EQ.2) kf=23    
        mstu10=mstu(10) 
        mstu(10)=1  
        p(nc,5)=ecmc    
        CALL lu1ent(nc,kf,ecmc,0.,0.)   
        k(nc,1)=21  
        k(nc,3)=1   
        mstu(10)=mstu10 
      ENDIF 
    
C...Choice of flavour and jet configuration.    
      CALL luxkfl(kfl,ecm,ecmc,kflc)    
      IF(kflc.EQ.0) goto 100    
      CALL luxjet(ecmc,njet,cut)    
      kfln=21   
      IF(njet.EQ.4) CALL lux4jt(njet,cut,kflc,ecmc,kfln,x1,x2,x4,   
     &x12,x14)  
      IF(njet.EQ.3) CALL lux3jt(njet,cut,kflc,ecmc,x1,x3)   
      IF(njet.EQ.2) mstj(120)=1 
    
C...Fill jet configuration and origin.  
      IF(njet.EQ.2.AND.mstj(101).NE.5) CALL lu2ent(nc+1,kflc,-kflc,ecmc)    
      IF(njet.EQ.2.AND.mstj(101).EQ.5) CALL lu2ent(-(nc+1),kflc,-kflc,  
     &ecmc) 
      IF(njet.EQ.3) CALL lu3ent(nc+1,kflc,21,-kflc,ecmc,x1,x3)  
      IF(njet.EQ.4.AND.kfln.EQ.21) CALL lu4ent(nc+1,kflc,kfln,kfln, 
     &-kflc,ecmc,x1,x2,x4,x12,x14)  
      IF(njet.EQ.4.AND.kfln.NE.21) CALL lu4ent(nc+1,kflc,-kfln,kfln,    
     &-kflc,ecmc,x1,x2,x4,x12,x14)  
      DO 110 ip=nc+1,n  
  110 k(ip,3)=k(ip,3)+min(mstj(115)/2,1)+(mstj(115)/3)*(nc-1)   
    
C...Angular orientation according to matrix element.    
      IF(mstj(106).EQ.1) THEN   
        CALL luxdif(nc,njet,kflc,ecmc,chi,the,phi)  
        CALL ludbrb(nc+1,n,0.,chi,0d0,0d0,0d0)  
        CALL ludbrb(nc+1,n,the,phi,0d0,0d0,0d0) 
      ENDIF 
    
C...Rotation and boost from radiative photon.   
      IF(mk.EQ.1) THEN  
        dbek=-pak/(ecm-pak) 
        nmin=nc+1-mstj(115)/3   
        CALL ludbrb(nmin,n,0.,-phik,0d0,0d0,0d0)    
        CALL ludbrb(nmin,n,alpk,0.,dbek*sin(thek),0d0,dbek*cos(thek))   
        CALL ludbrb(nmin,n,0.,phik,0d0,0d0,0d0) 
      ENDIF 
    
C...Generate parton shower. Rearrange along strings and check.  
      IF(mstj(101).EQ.5) THEN   
        CALL lushow(n-1,n,ecmc) 
        mstj14=mstj(14) 
        IF(mstj(105).EQ.-1) mstj(14)=0  
        IF(mstj(105).GE.0) mstu(28)=0   
        CALL luprep(0)  
        mstj(14)=mstj14 
        IF(mstj(105).GE.0.AND.mstu(28).NE.0) goto 100   
      ENDIF 
    
C...Fragmentation/decay generation. Information for LUTABU. 
      IF(mstj(105).EQ.1) CALL luexec    
      mstu(161)=kflc    
      mstu(162)=-kflc   
    
      RETURN    
      END