Update from Jiri: enables the task to create a very simple offline L1 jet trigger...

[u/mrichter/AliRoot.git] / HIJING / hijing1_36 / hijing.F

diff --git a/HIJING/hijing1_36/hijing.F b/HIJING/hijing1_36/hijing.F
index 89ad5c9a7408259285d87415491ef9bbcda17ea9..07db9ceb025dd87e8454cdf3b87c414cd853e03f 100644 (file)
--- a/HIJING/hijing1_36/hijing.F
+++ b/HIJING/hijing1_36/hijing.F
@@ -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,16 +377,14 @@ 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
 
-#if 1
 c *** cl glauber ***
-       IF(NCOLT.EQ.0) THEN 
-           goto 60
-        end if
-
         npart=0
         xmeana=0D0
         ymeana=0D0
@@ -397,72 +396,74 @@ c *** cl glauber ***
         ym2=0D0
         xym=0D0
 
-       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(SCIP(JP,JT).GT.-1.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(SCIP(JP,JT).GT.-1.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(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
 
-        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
+           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
 
-        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,' ',eccrp,
-     1  ' ',eccpart
-        RETURN
-#endif
+           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.