]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - HIJING/hijing1_36/hijing.F
Update from Jiri: enables the task to create a very simple offline L1 jet trigger...
[u/mrichter/AliRoot.git] / HIJING / hijing1_36 / hijing.F
index 9900358a722f59be5a5761df61898642128b6da0..07db9ceb025dd87e8454cdf3b87c414cd853e03f 100644 (file)
@@ -154,7 +154,7 @@ C****************************************************************
        SUBROUTINE HIJING(FRAME,BMIN0,BMAX0)
        CHARACTER FRAME*8
        DIMENSION SCIP(300,300),RNIP(300,300),SJIP(300,300),JTP(3),
-     &                 IPCOL(90000),ITCOL(90000)
+     &                 IPCOL(90000),ITCOL(90000),SCIP2(300,300)
 #define BLANKET_SAVE
 #include "hiparnt.inc"
 C
@@ -342,6 +342,7 @@ C
           SCIP(JP,JT)=-1.0
           B2=(YP(1,JP)+BBX-YT(1,JT))**2+(YP(2,JP)+BBY-YT(2,JT))**2
           R2=B2*HIPR1(40)/HIPR1(31)/0.1
+          SCIP2(JP,JT)=R2
 C              ********mb=0.1*fm, YP is in fm,HIPR1(31) is in mb
           RRB1=MIN((YP(1,JP)**2+YP(2,JP)**2)
      &          /1.2**2/REAL(IHNT2(1))**0.6666667,1.0)
@@ -376,11 +377,93 @@ C                 ********perform elastic collisions
           RNIP(JP,JT)=RANTOT
           SJIP(JP,JT)=HINT1(18)
           NCOLT=NCOLT+1
+           if (R2.GT.2.D0) THEN
+              write (8,*) R2
+           ENDIF
           IPCOL(NCOLT)=JP
           ITCOL(NCOLT)=JT
 70     CONTINUE
+
+c *** cl glauber ***
+        npart=0
+        xmeana=0D0
+        ymeana=0D0
+        xmeanb=0D0
+        ymeanb=0D0
+        xmeanp=0D0
+        ymeanp=0D0
+        xm2=0D0
+        ym2=0D0
+        xym=0D0
+
+       IF(NCOLT>0) THEN 
+           DO 1110 JP=1,IHNT2(1)
+              YP(1,JP)=YP(1,JP)+BBX
+              YP(2,JP)=YP(2,JP)+BBY
+              xmeana=xmeana+YP(1,JP)
+              ymeana=ymeana+YP(2,JP)
+              DO 1120 JT=1,IHNT2(3)
+                 IF(SCIP2(JP,JT).LT.2.0D0) THEN
+                    npart=npart+1
+                    xmeanp=xmeanp+YP(1,JP)
+                    ymeanp=ymeanp+YP(2,JP)
+                    xm2=xm2+YP(1,JP)*YP(1,JP)
+                    ym2=ym2+YP(2,JP)*YP(2,JP)
+                    xym=xym+YP(1,JP)*YP(2,JP)
+                    goto 1110
+                 end if
+ 1120         continue
+ 1110      continue
+
+           DO 1130 JT=1,IHNT2(3)
+              xmeanb=xmeanb+YT(1,JT)
+              ymeanb=ymeanb+YT(2,JT)
+              DO 1140 JP=1,IHNT2(1)
+                 IF(SCIP2(JP,JT).LT.2.0D0) THEN
+                    npart=npart+1
+                    xmeanp=xmeanp+YT(1,JT)
+                    ymeanp=ymeanp+YT(2,JT)
+                    xm2=xm2+YT(1,JT)*YT(1,JT)
+                    ym2=ym2+YT(2,JT)*YT(2,JT)
+                    xym=xym+YT(1,JT)*YT(2,JT)
+                    goto 1130
+                 end if
+ 1140         continue
+ 1130      continue
+
+           IF (npart.GT.0) THEN
+              xmeana=xmeana/IHNT2(1)
+              ymeana=ymeana/IHNT2(1)
+              xmeanb=xmeanb/IHNT2(3)
+              ymeanb=ymeanb/IHNT2(3)
+              xmeanp=xmeanp/npart
+              ymeanp=ymeanp/npart
+              xm2=xm2/npart
+              ym2=ym2/npart
+              xym=xym/npart
+
+              sx2=xm2-xmeanp*xmeanp
+              sy2=ym2-ymeanp*ymeanp
+              sxy=xym-xmeanp*ymeanp
+           
+              delx=xmeanb-xmeana
+              dely=ymeanb-ymeana
+              dtmp=delx**2+dely**2
+              bbtrue=sqrt(dtmp)
+              dnumt=(sy2-sx2)*(delx**2-dely**2)-4D0*sxy*delx*dely
+              ddent=(sy2+sx2)*bbtrue**2
+              eccrp=dnumt/ddent
+              dtmp=(sy2-sx2)*(sy2-sx2)+4D0*sxy*sxy
+              eccpart=sqrt(dtmp)/(sx2+sy2)
+              eccmc=(sy2-sx2)/(sy2+sx2)
+              write(*,*),'HOUT: ',bb,' ',bbtrue,' ',ncolt,' ',npart,
+     1             ' ',eccrp,' ',eccpart
+           end if
+        end if
+
 C              ********total number interactions proj and targ has
 C                              suffered
+
        IF(NCOLT.EQ.0) THEN
           NLOP=NLOP+1
            IF(NLOP.LE.20.OR.
@@ -601,8 +684,6 @@ 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
@@ -623,10 +704,10 @@ C
                IF(IHPR2(21).EQ.0) THEN
                   CALL LUEDIT_HIJING(2)
                ELSE
-351               N_ST=N_ST+1
+351                N_ST=N_ST+1
                   IF(K(N_ST,2).LT.91.OR.K(N_ST,2).GT.93) GO TO  351
-                  IDSTR=K(N_ST,2)
-                  N_ST=N_ST+1
+                   IDSTR=K(N_ST,2)
+                   N_ST=N_ST+1
                ENDIF
 C
                IF(FRAME.EQ.'LAB') THEN
@@ -659,31 +740,12 @@ 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