lu2ent.f

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

    
C*********************************************************************  
    
      SUBROUTINE lu2ent(IP,KF1,KF2,PECM)    
    
C...Purpose: to store two partons/particles in their CM frame,  
C...with the first along the +z axis.   
      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...Standard checks.    
      mstu(28)=0    
      IF(mstu(12).GE.1) CALL lulist(0)  
      ipa=max(1,iabs(ip))   
      IF(ipa.GT.mstu(4)-1) CALL luerrm(21,  
     &'(LU2ENT:) writing outside LUJETS memory')    
      kc1=lucomp(kf1)   
      kc2=lucomp(kf2)   
      IF(kc1.EQ.0.OR.kc2.EQ.0) CALL luerrm(12,  
     &'(LU2ENT:) unknown flavour code') 
    
C...Find masses. Reset K, P and V vectors.  
      pm1=0.    
      IF(mstu(10).EQ.1) pm1=p(ipa,5)    
      IF(mstu(10).GE.2) pm1=ulmass(kf1) 
      pm2=0.    
      IF(mstu(10).EQ.1) pm2=p(ipa+1,5)  
      IF(mstu(10).GE.2) pm2=ulmass(kf2) 
      DO 100 i=ipa,ipa+1    
      DO 100 j=1,5  
      k(i,j)=0  
      p(i,j)=0. 
  100 v(i,j)=0. 
    
C...Check flavours. 
      kq1=kchg(kc1,2)*isign(1,kf1)  
      kq2=kchg(kc2,2)*isign(1,kf2)  
      IF(kq1+kq2.NE.0.AND.kq1+kq2.NE.4) CALL luerrm(2,  
     &'(LU2ENT:) unphysical flavour combination')   
      k(ipa,2)=kf1  
      k(ipa+1,2)=kf2    
    
C...Store partons/particles in K vectors for normal case.   
      IF(ip.GE.0) THEN  
        k(ipa,1)=1  
        IF(kq1.NE.0.AND.kq2.NE.0) k(ipa,1)=2    
        k(ipa+1,1)=1    
    
C...Store partons in K vectors for parton shower evolution. 
      ELSE  
        IF(kq1.EQ.0.OR.kq2.EQ.0) CALL luerrm(2, 
     &  '(LU2ENT:) requested flavours can not develop parton shower')   
        k(ipa,1)=3  
        k(ipa+1,1)=3    
        k(ipa,4)=mstu(5)*(ipa+1)    
        k(ipa,5)=k(ipa,4)   
        k(ipa+1,4)=mstu(5)*ipa  
        k(ipa+1,5)=k(ipa+1,4)   
      ENDIF 
    
C...Check kinematics and store partons/particles in P vectors.  
      IF(pecm.LE.pm1+pm2) CALL luerrm(13,   
     &'(LU2ENT:) energy smaller than sum of masses')    
      pa=sqrt(max(0.,(pecm**2-pm1**2-pm2**2)**2-(2.*pm1*pm2)**2))/  
     &(2.*pecm) 
      p(ipa,3)=pa   
      p(ipa,4)=sqrt(pm1**2+pa**2)   
      p(ipa,5)=pm1  
      p(ipa+1,3)=-pa    
      p(ipa+1,4)=sqrt(pm2**2+pa**2) 
      p(ipa+1,5)=pm2    
    
C...Set N. Optionally fragment/decay.   
      n=ipa+1   
      IF(ip.EQ.0) CALL luexec   
    
      RETURN    
      END