attrad.f

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

C
C
C
C****************************************************************
C       conduct soft radiation according to dipole approxiamtion
C****************************************************************
        SUBROUTINE attrad(IERROR)
C
        common/hiparnt/hipr1(100),ihpr2(50),hint1(100),ihnt2(50)
        SAVE  /hiparnt/
        common/hijdat/hidat0(10,10),hidat(10)
        SAVE  /hijdat/
        common/lujets/n,k(9000,5),p(9000,5),v(9000,5)
        SAVE  /lujets/
        ierror=0

C.....S INVARIANT MASS-SQUARED BETWEEN PARTONS I AND I+1......
C.....SM IS THE LARGEST MASS-SQUARED....

40      sm=0.
        jl=1
        DO 30 i=1,n-1
           s=2.*(p(i,4)*p(i+1,4)-p(i,1)*p(i+1,1)-p(i,2)*p(i+1,2)
     &          -p(i,3)*p(i+1,3))+p(i,5)**2+p(i+1,5)**2
           IF(s.LT.0.) s=0.
           wp=sqrt(s)-1.5*(p(i,5)+p(i+1,5))
           IF(wp.GT.sm) THEN
              pbt1=p(i,1)+p(i+1,1)
              pbt2=p(i,2)+p(i+1,2)
              pbt3=p(i,3)+p(i+1,3)
              pbt4=p(i,4)+p(i+1,4)
              btt=(pbt1**2+pbt2**2+pbt3**2)/pbt4**2
              IF(btt.GE.1.0-1.0e-10) go to 30
              IF((i.NE.1.OR.i.NE.n-1).AND.
     &             (k(i,2).NE.21.AND.k(i+1,2).NE.21)) go to 30
              jl=i
              sm=wp
           ENDIF
30      CONTINUE
        s=(sm+1.5*(p(jl,5)+p(jl+1,5)))**2
        IF(sm.LT.hipr1(5)) goto 2
     
C.....MAKE PLACE FOR ONE GLUON.....
        IF(jl+1.EQ.n) goto 190
        DO 160 j=n,jl+2,-1
                k(j+1,1)=k(j,1)
                k(j+1,2)=k(j,2)
                DO 150 m=1,5
C+++BAC
                   v(j+1,m) = v(j,m)
C---BAC
150                     p(j+1,m)=p(j,m)
160             CONTINUE
190     n=n+1
     
C.....BOOST TO REST SYSTEM FOR PARTICLES JL AND JL+1.....
        p1=p(jl,1)+p(jl+1,1)
        p2=p(jl,2)+p(jl+1,2)
        p3=p(jl,3)+p(jl+1,3)
        p4=p(jl,4)+p(jl+1,4)
        bex=-p1/p4
        bey=-p2/p4
        bez=-p3/p4
        imin=jl
        imax=jl+1
        CALL atrobo(0.,0.,bex,bey,bez,imin,imax,ierror)
        IF(ierror.NE.0) RETURN
C.....ROTATE TO Z-AXIS....
        cth=p(jl,3)/sqrt(p(jl,4)**2-p(jl,5)**2)
        IF(abs(cth).GT.1.0)  cth=max(-1.,min(1.,cth))
        theta=acos(cth)
        phi=ulangl(p(jl,1),p(jl,2))
        CALL atrobo(0.,-phi,0.,0.,0.,imin,imax,ierror)
        CALL atrobo(-theta,0.,0.,0.,0.,imin,imax,ierror)
     
C.....CREATE ONE GLUON AND ORIENTATE.....
     
1       CALL ar3jet(s,x1,x3,jl)
        CALL arorie(s,x1,x3,jl)         
        IF(hidat(2).GT.0.0) THEN
           ptg1=sqrt(p(jl,1)**2+p(jl,2)**2)
           ptg2=sqrt(p(jl+1,1)**2+p(jl+1,2)**2)
           ptg3=sqrt(p(jl+2,1)**2+p(jl+2,2)**2)
           ptg=max(ptg1,ptg2,ptg3)
           IF(ptg.GT.hidat(2)) THEN
              fmfact=exp(-(ptg**2-hidat(2)**2)/hipr1(2)**2)
              IF(atl_ran(nseed).GT.fmfact) go to 1
           ENDIF
        ENDIF
C.....ROTATE AND BOOST BACK.....
        imin=jl
        imax=jl+2
        CALL atrobo(theta,phi,-bex,-bey,-bez,imin,imax,ierror)
        IF(ierror.NE.0) RETURN
C.....ENUMERATE THE GLUONS.....
        k(jl+2,1)=k(jl+1,1)
        k(jl+2,2)=k(jl+1,2)
        k(jl+2,3)=k(jl+1,3)
        k(jl+2,4)=k(jl+1,4)
        k(jl+2,5)=k(jl+1,5)
        p(jl+2,5)=p(jl+1,5)
        k(jl+1,1)=2
        k(jl+1,2)=21
        k(jl+1,3)=0
        k(jl+1,4)=0
        k(jl+1,5)=0
        p(jl+1,5)=0.

        v(jl+1,1) = 0.
        v(jl+1,2) = 0.
        v(jl+1,3) = 0.
        v(jl+1,4) = 0.
        v(jl+1,5) = 0.

C----THETA FUNCTION DAMPING OF THE EMITTED GLUONS. FOR HADRON-HADRON.
C----R0=VFR(2)
C       IF(VFR(2).GT.0.) THEN
C       PTG=SQRT(P(JL+1,1)**2+P(JL+1,2)**2)
C       PTGMAX=WSTRI/2.
C       DOPT=SQRT((4.*PAR(71)*VFR(2))/WSTRI)
C       PTOPT=(DOPT*WSTRI)/(2.*VFR(2))
C       IF(PTG.GT.PTOPT) IORDER=IORDER-1
C       IF(PTG.GT.PTOPT) GOTO 1
C       ENDIF
C-----
        IF(sm.GE.hipr1(5)) goto 40

2       k(1,1)=2
        k(1,3)=0
        k(1,4)=0
        k(1,5)=0
        k(n,1)=1
        k(n,3)=0
        k(n,4)=0
        k(n,5)=0

        RETURN
        END