]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - HIJING/hijing1_36/hijing.F
Correct ran function.
[u/mrichter/AliRoot.git] / HIJING / hijing1_36 / hijing.F
index fb077dcbcc4cbfd052030306e36fa32c9b68535b..9df2e9215acce7fd7722fba4a6c5f7e9ddf33a27 100644 (file)
@@ -159,6 +159,7 @@ C****************************************************************
 #include "hiparnt.inc"
 C
 #include "hijcrdn.inc"
+#include "hijglbr.inc"
 #include "himain1.inc"
 #include "himain2.inc"
 #include "histrng.inc"
@@ -186,6 +187,15 @@ C
        IF(IHNT2(1).LE.1) GO TO 14
        DO 10 KP=1,IHNT2(1)
 5      R=HIRND(1)
+c
+        if(IHNT2(1).EQ.2) then
+           rnd1=max(RLU_HIJING(NSEED),1.0e-20)
+           rnd2=max(RLU_HIJING(NSEED),1.0e-20)
+           rnd3=max(RLU_HIJING(NSEED),1.0e-20)
+           R=-0.5*(log(rnd1)*4.38/2.0+log(rnd2)*0.85/2.0
+     &          +4.38*0.85*log(rnd3)/(4.38+0.85))
+        endif
+c
        X=RLU_HIJING(0)
        CX=2.0*X-1.0
        SX=SQRT(1.0-CX*CX)
@@ -206,6 +216,13 @@ C                  ********two neighbors cannot be closer than
 C                              HIPR1(29)
 8      CONTINUE
 10     CONTINUE
+c*******************************
+        if(IHNT2(1).EQ.2) then
+           YP(1,2)=-YP(1,1)
+           YP(2,2)=-YP(2,1)
+           YP(3,2)=-YP(3,1)
+        endif
+c********************************
        DO 12 I=1,IHNT2(1)-1
        DO 12 J=I+1,IHNT2(1)
        IF(YP(3,I).GT.YP(3,J)) GO TO 12
@@ -227,6 +244,15 @@ C******************************
        IF(IHNT2(3).LE.1) GO TO 24
        DO 20 KT=1,IHNT2(3)
 15     R=HIRND(2)
+c
+         if(IHNT2(3).EQ.2) then
+            rnd1=max(RLU_HIJING(NSEED),1.0e-20)
+            rnd2=max(RLU_HIJING(NSEED),1.0e-20)
+            rnd3=max(RLU_HIJING(NSEED),1.0e-20)
+            R=-0.5*(log(rnd1)*4.38/2.0+log(rnd2)*0.85/2.0
+     &           +4.38*0.85*log(rnd3)/(4.38+0.85))
+         endif
+c
        X=RLU_HIJING(0)
        CX=2.0*X-1.0
        SX=SQRT(1.0-CX*CX)
@@ -247,6 +273,13 @@ C                  ********two neighbors cannot be closer than
 C                              HIPR1(29)
 18     CONTINUE
 20     CONTINUE
+c**********************************
+         if(IHNT2(3).EQ.2) then
+            YT(1,2)=-YT(1,1)
+            YT(2,2)=-YT(2,1)
+            YT(3,2)=-YT(3,1)
+         endif
+c*********************************
        DO 22 I=1,IHNT2(3)-1
        DO 22 J=I+1,IHNT2(3)
        IF(YT(3,I).LT.YT(3,J)) GO TO 22
@@ -282,6 +315,10 @@ C                  ********Initialize for a new event
        N01=0
        N10=0
        N11=0
+        NELT=0
+        NINT=0
+        NELP=0
+        NINP=0
        NSG=0
        NCOLT=0
 
@@ -501,6 +538,26 @@ C          ********conduct soft scattering between JP and JT
 
 200    CONTINUE
 
+c
+c**************************
+c
+       DO 201 JP=1,IHNT2(1)
+           IF(NFP(JP,5).GT.2) THEN
+              NINP=NINP+1
+           ELSE IF(NFP(JP,5).EQ.2.OR.NFP(JP,5).EQ.1) THEN
+              NELP=NELP+1
+           ENDIF
+ 201    continue
+       DO 202 JT=1,IHNT2(3)
+           IF(NFT(JT,5).GT.2) THEN
+              NINT=NINT+1
+           ELSE IF(NFT(JT,5).EQ.2.OR.NFT(JT,5).EQ.1) THEN
+              NELT=NELT+1
+           ENDIF
+ 202    continue
+c     
+c*******************************
+
 C********perform jet quenching for jets with PT>HIPR1(11)**********
 
        IF((IHPR2(8).NE.0.OR.IHPR2(3).NE.0).AND.IHPR2(4).GT.0.AND.
@@ -523,6 +580,9 @@ C********N_STR+1 is the number of strings in fragmentation
 C********the number of strings before a line is stored in K(I,4)
 C********IDSTR is id number of the string system (91,92 or 93)
 C
+        GAMMA0=HINT1(1)/2.0/0.93827
+        BETA0=sqrt(GAMMA0**2-1.0)/GAMMA0
+c
         IF(IHPR2(20).NE.0) THEN
           DO 360 ISG=1,NSG
                CALL HIJFRG(ISG,3,IERROR)
@@ -574,12 +634,31 @@ C       ****** identify the mother particle
                   PATT(NATT,2)=P(I,2)
                   PATT(NATT,3)=P(I,3)
                   PATT(NATT,4)=P(I,4)
-                   VATT(NATT,1)=V(I,1)
-                   VATT(NATT,2)=V(I,2)
-                   VATT(NATT,3)=V(I,3)
-                   VATT(NATT,4)=V(I,4)
-                   
                   EATT=EATT+P(I,4)
+                   VATT01=0.5*(YP(1,IASG(ISG,1))+YT(1,IASG(ISG,2)))
+                   VATT02=0.5*(YP(2,IASG(ISG,1))+YT(2,IASG(ISG,2)))
+                   VATT03=0.5*(YP(3,IASG(ISG,1))+YT(3,IASG(ISG,2)))
+     &                  /GAMMA0
+                   VATT04=-0.5*(YP(3,IASG(ISG,1))
+     &                  -YT(3,IASG(ISG,2)))/GAMMA0/BETA0
+                   RARB=1.12*(IHNT2(1)**0.33333+IHNT2(3)**0.33333)
+                   V3MIN1=RARB/GAMMA0
+                   V3MIN2=1.0/MAX(1.0,5.08*ABS(PATT(I,3)))
+                   VATT_MIN=MAX(V3MIN1,V3MIN2)
+                   VATT03=VATT03+(0.5-RLU_HIJING(0))*VATT_MIN
+                   amt2=P(I,1)**2+P(I,2)**2+P(I,5)**2
+                   IF(amt2.GT.0.0) THEN
+                      tauf=0.2*2.0*P(I,3)/amt2
+                      VATT(NATT,1)=VATT01
+                      VATT(NATT,2)=VATT02
+                      VATT(NATT,3)=VATT03+tauf
+                      VATT(NATT,4)=abs(VATT(NATT,3))
+                   ELSE
+                      VATT(NATT,4)=abs(VATT03)
+                      VATT(NATT,1)=VATT01
+                      VATT(NATT,2)=VATT02
+                      VATT(NATT,3)=VATT03
+                   ENDIF
 360       CONTINUE
 C              ********Fragment the q-qbar jets systems *****
 C
@@ -602,12 +681,9 @@ C
                N_ST=1
                IDSTR=92
 
-               NFTP=NFP(J_JTP,5)
-               IF(NTP.EQ.2) NFTP=10+NFT(J_JTP,5)
-
                IF(IHPR2(21).EQ.0) THEN
                   CALL LUEDIT_HIJING(2)
-               ELSE IF (NFTP.EQ. 3 .OR. NFTP .EQ. 13) THEN
+               ELSE
 381               N_ST=N_ST+1
                   IF(K(N_ST,2).LT.91.OR.K(N_ST,2).GT.93) GO TO  381
                   IDSTR=K(N_ST,2)
@@ -618,6 +694,8 @@ C
                ENDIF
 C              ******** boost back to lab frame(if it was in)
 C
+               NFTP=NFP(J_JTP,5)
+               IF(NTP.EQ.2) NFTP=10+NFT(J_JTP,5)
                N_STR=0
                DO 390 I=N_ST,N
                   IF(K(I,2).EQ.IDSTR) THEN
@@ -629,26 +707,46 @@ C
                   KATT(NATT,1)=K(I,2)
                   KATT(NATT,2)=NFTP
                   KATT(NATT,4)=K(I,1)
-                  IF(K(I,3).EQ.0 ) THEN
+                   IF(K(I,3).EQ.0 .OR. K(K(I,3),2).EQ.IDSTR) THEN
                      KATT(NATT,3)=0
-                   ELSE
-                      IF(K(K(I,3),2).EQ.IDSTR) THEN
-                         KATT(NATT,3)=0
-                      ELSE
-                         KATT(NATT,3)=NATT-I+K(I,3)+N_STR-K(K(I,3),4)
-                      ENDIF
+                  ELSE
+                     KATT(NATT,3)=NATT-I+K(I,3)+N_STR-K(K(I,3),4)
                   ENDIF
 C       ****** identify the mother particle
                   PATT(NATT,1)=P(I,1)
                   PATT(NATT,2)=P(I,2)
                   PATT(NATT,3)=P(I,3)
                   PATT(NATT,4)=P(I,4)
-                   VATT(NATT,1)=V(I,1)
-                   VATT(NATT,2)=V(I,2)
-                   VATT(NATT,3)=V(I,3)
-                   VATT(NATT,4)=V(I,4)
-
                   EATT=EATT+P(I,4)
+C     ******* the following calculate the formation time and position
+C             of the produced hadrons
+                   IF(NTP.EQ.1) THEN
+                      VATT01=YP(1,J_JTP)
+                      VATT02=YP(2,J_JTP)
+                      VATT03=YP(3,J_JTP)/GAMMA0
+                   ELSE IF(NTP.EQ.2) THEN
+                      VATT01=YT(1,J_JTP)
+                      VATT02=YT(2,J_JTP)
+                      VATT03=YT(3,J_JTP)/GAMMA0
+                   ENDIF
+                   RARB=1.12*(IHNT2(1)**0.33333+IHNT2(3)**0.33333)
+                   V3MIN1=RARB/GAMMA0
+                   V3MIN2=1.0/MAX(1.0,5.08*ABS(PATT(I,3)))
+                   VATT_MIN=MAX(V3MIN1,V3MIN2)
+                   VATT03=VATT03+(0.5-RLU_HIJING(0))*VATT_MIN
+                   amt2=P(I,1)**2+P(I,2)**2+P(I,5)**2
+                   IF(amt2.GT.0.0) THEN
+                      tauf=0.2*2.0*P(I,3)/amt2
+                      VATT(NATT,1)=VATT01
+                      VATT(NATT,2)=VATT02
+                      VATT(NATT,3)=VATT03+tauf
+                      VATT(NATT,4)=abs(VATT(NATT,3))
+                   ELSE
+                      VATT(NATT,4)=abs(VATT03)
+                      VATT(NATT,1)=VATT01
+                      VATT(NATT,2)=VATT02
+                      VATT(NATT,3)=VATT03
+                   ENDIF
 390            CONTINUE 
 400       CONTINUE
 C              ********Fragment the q-qq related string systems
@@ -663,12 +761,20 @@ C         ********Fragment the q-qq related string systems
                PATT(NATT,2)=PDR(I,2)
                PATT(NATT,3)=PDR(I,3)
                PATT(NATT,4)=PDR(I,4)
-                VATT(NATT,1)=V(I,1)
-                VATT(NATT,2)=V(I,2)
-                VATT(NATT,3)=V(I,3)
-                VATT(NATT,4)=V(I,4)
-
                EATT=EATT+PDR(I,4)
+                amt2=P(I,1)**2+P(I,2)**2+P(I,5)**2
+                IF(amt2.GT.0.0) THEN
+                   tauf=0.2*2.0*P(I,3)/amt2
+                   VATT(NATT,1)=0.0
+                   VATT(NATT,2)=0.0
+                   VATT(NATT,3)=tauf
+                   VATT(NATT,4)=abs(VATT(NATT,3))
+                ELSE
+                   VATT(NATT,4)=0.0
+                   VATT(NATT,1)=0.0
+                   VATT(NATT,2)=0.0
+                   VATT(NATT,3)=0.0
+                ENDIF
 450    CONTINUE
 C                      ********store the direct-produced particles
 C