author morsch Wed, 12 Sep 2012 12:57:56 +0000 (12:57 +0000) committer morsch Wed, 12 Sep 2012 12:57:56 +0000 (12:57 +0000)

@@ -380,12 +380,7 @@ C                  ********perform elastic collisions
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,69 +392,68 @@ 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(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

-        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(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

-        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
+           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

C              ********total number interactions proj and targ has
C                              suffered