HIJING/hipyset1_35/lucell_hijing.F

   1 * $Id$
   2
   3 C*********************************************************************
   4
   5       SUBROUTINE LUCELL_HIJING(NJET)
   6
   7 C...Purpose: to provide a simple way of jet finding in an eta-phi-ET
   8 C...coordinate frame, as used for calorimeters at hadron colliders.
   9 #include "lujets_hijing.inc"
  10 #include "ludat1_hijing.inc"
  11 #include "ludat2_hijing.inc"
  12
  13 C...Loop over all particles. Find cell that was hit by given particle.
  14       NCE2=2*MSTU(51)*MSTU(52)
  15       PTLRAT=1./SINH(PARU(51))**2
  16       NP=0
  17       NC=N
  18       DO 110 I=1,N
  19       IF(K(I,1).LE.0.OR.K(I,1).GT.10) GOTO 110
  20       IF(P(I,1)**2+P(I,2)**2.LE.PTLRAT*P(I,3)**2) GOTO 110
  21       IF(MSTU(41).GE.2) THEN
  22         KC=LUCOMP_HIJING(K(I,2))
  23         IF(KC.EQ.0.OR.KC.EQ.12.OR.KC.EQ.14.OR.KC.EQ.16.OR.
  24      &  KC.EQ.18) GOTO 110
  25         IF(MSTU(41).GE.3.AND.KCHG(KC,2).EQ.0.AND.LUCHGE_HIJING(K(I,2))
  26      $       .EQ.0)GOTO 110
  27       ENDIF
  28       NP=NP+1
  29       PT=SQRT(P(I,1)**2+P(I,2)**2)
  30       ETA=SIGN(LOG((SQRT(PT**2+P(I,3)**2)+ABS(P(I,3)))/PT),P(I,3))
  31       IETA=MAX(1,MIN(MSTU(51),1+INT(MSTU(51)*0.5*(ETA/PARU(51)+1.))))
  32       PHI=ULANGL_HIJING(P(I,1),P(I,2))
  33       IPHI=MAX(1,MIN(MSTU(52),1+INT(MSTU(52)*0.5*(PHI/PARU(1)+1.))))
  34       IETPH=MSTU(52)*IETA+IPHI
  35
  36 C...Add to cell already hit, or book new cell.
  37       DO 100 IC=N+1,NC
  38       IF(IETPH.EQ.K(IC,3)) THEN
  39         K(IC,4)=K(IC,4)+1
  40         P(IC,5)=P(IC,5)+PT
  41         GOTO 110
  42       ENDIF
  43   100 CONTINUE
  44       IF(NC.GE.MSTU(4)-MSTU(32)-5) THEN
  45          CALL LUERRM_HIJING(11
  46      $        ,'(LUCELL_HIJING:) no more memory left in LUJETS_HIJING')
  47         NJET=-2
  48         RETURN
  49       ENDIF
  50       NC=NC+1
  51       K(NC,3)=IETPH
  52       K(NC,4)=1
  53       K(NC,5)=2
  54       P(NC,1)=(PARU(51)/MSTU(51))*(2*IETA-1-MSTU(51))
  55       P(NC,2)=(PARU(1)/MSTU(52))*(2*IPHI-1-MSTU(52))
  56       P(NC,5)=PT
  57   110 CONTINUE
  58
  59 C...Smear true bin content by calorimeter resolution.
  60       IF(MSTU(53).GE.1) THEN
  61         DO 130 IC=N+1,NC
  62         PEI=P(IC,5)
  63         IF(MSTU(53).EQ.2) PEI=P(IC,5)/COSH(P(IC,1))
  64   120   PEF=PEI+PARU(55)*SQRT(-2.*LOG(MAX(1E-10,RLU_HIJING(0)))*PEI)*
  65      &  COS(PARU(2)*RLU_HIJING(0))
  66         IF(PEF.LT.0..OR.PEF.GT.PARU(56)*PEI) GOTO 120
  67         P(IC,5)=PEF
  68   130   IF(MSTU(53).EQ.2) P(IC,5)=PEF*COSH(P(IC,1))
  69       ENDIF
  70
  71 C...Find initiator cell: the one with highest pT of not yet used ones.
  72       NJ=NC
  73   140 ETMAX=0.
  74       DO 150 IC=N+1,NC
  75       IF(K(IC,5).NE.2) GOTO 150
  76       IF(P(IC,5).LE.ETMAX) GOTO 150
  77       ICMAX=IC
  78       ETA=P(IC,1)
  79       PHI=P(IC,2)
  80       ETMAX=P(IC,5)
  81   150 CONTINUE
  82       IF(ETMAX.LT.PARU(52)) GOTO 210
  83       IF(NJ.GE.MSTU(4)-MSTU(32)-5) THEN
  84          CALL LUERRM_HIJING(11
  85      $        ,'(LUCELL_HIJING:) no more memory left in LUJETS_HIJING')
  86         NJET=-2
  87         RETURN
  88       ENDIF
  89       K(ICMAX,5)=1
  90       NJ=NJ+1
  91       K(NJ,4)=0
  92       K(NJ,5)=1
  93       P(NJ,1)=ETA
  94       P(NJ,2)=PHI
  95       P(NJ,3)=0.
  96       P(NJ,4)=0.
  97       P(NJ,5)=0.
  98
  99 C...Sum up unused cells within required distance of initiator.
 100       DO 160 IC=N+1,NC
 101       IF(K(IC,5).EQ.0) GOTO 160
 102       IF(ABS(P(IC,1)-ETA).GT.PARU(54)) GOTO 160
 103       DPHIA=ABS(P(IC,2)-PHI)
 104       IF(DPHIA.GT.PARU(54).AND.DPHIA.LT.PARU(2)-PARU(54)) GOTO 160
 105       PHIC=P(IC,2)
 106       IF(DPHIA.GT.PARU(1)) PHIC=PHIC+SIGN(PARU(2),PHI)
 107       IF((P(IC,1)-ETA)**2+(PHIC-PHI)**2.GT.PARU(54)**2) GOTO 160
 108       K(IC,5)=-K(IC,5)
 109       K(NJ,4)=K(NJ,4)+K(IC,4)
 110       P(NJ,3)=P(NJ,3)+P(IC,5)*P(IC,1)
 111       P(NJ,4)=P(NJ,4)+P(IC,5)*PHIC
 112       P(NJ,5)=P(NJ,5)+P(IC,5)
 113   160 CONTINUE
 114
 115 C...Reject cluster below minimum ET, else accept.
 116       IF(P(NJ,5).LT.PARU(53)) THEN
 117         NJ=NJ-1
 118         DO 170 IC=N+1,NC
 119   170   IF(K(IC,5).LT.0) K(IC,5)=-K(IC,5)
 120       ELSEIF(MSTU(54).LE.2) THEN
 121         P(NJ,3)=P(NJ,3)/P(NJ,5)
 122         P(NJ,4)=P(NJ,4)/P(NJ,5)
 123         IF(ABS(P(NJ,4)).GT.PARU(1)) P(NJ,4)=P(NJ,4)-SIGN(PARU(2),
 124      &  P(NJ,4))
 125         DO 180 IC=N+1,NC
 126   180   IF(K(IC,1).LT.0) K(IC,1)=0
 127       ELSE
 128         DO 190 J=1,4
 129   190   P(NJ,J)=0.
 130         DO 200 IC=N+1,NC
 131         IF(K(IC,5).GE.0) GOTO 200
 132         P(NJ,1)=P(NJ,1)+P(IC,5)*COS(P(IC,2))
 133         P(NJ,2)=P(NJ,2)+P(IC,5)*SIN(P(IC,2))
 134         P(NJ,3)=P(NJ,3)+P(IC,5)*SINH(P(IC,1))
 135         P(NJ,4)=P(NJ,4)+P(IC,5)*COSH(P(IC,1))
 136         K(IC,5)=0
 137   200   CONTINUE
 138       ENDIF
 139       GOTO 140
 140
 141 C...Arrange clusters in falling ET sequence.
 142   210 DO 230 I=1,NJ-NC
 143       ETMAX=0.
 144       DO 220 IJ=NC+1,NJ
 145       IF(K(IJ,5).EQ.0) GOTO 220
 146       IF(P(IJ,5).LT.ETMAX) GOTO 220
 147       IJMAX=IJ
 148       ETMAX=P(IJ,5)
 149   220 CONTINUE
 150       K(IJMAX,5)=0
 151       K(N+I,1)=31
 152       K(N+I,2)=98
 153       K(N+I,3)=I
 154       K(N+I,4)=K(IJMAX,4)
 155       K(N+I,5)=0
 156       DO 230 J=1,5
 157       P(N+I,J)=P(IJMAX,J)
 158   230 V(N+I,J)=0.
 159       NJET=NJ-NC
 160
 161 C...Convert to massless or massive four-vectors.
 162       IF(MSTU(54).EQ.2) THEN
 163         DO 240 I=N+1,N+NJET
 164         ETA=P(I,3)
 165         P(I,1)=P(I,5)*COS(P(I,4))
 166         P(I,2)=P(I,5)*SIN(P(I,4))
 167         P(I,3)=P(I,5)*SINH(ETA)
 168         P(I,4)=P(I,5)*COSH(ETA)
 169   240   P(I,5)=0.
 170       ELSEIF(MSTU(54).GE.3) THEN
 171         DO 250 I=N+1,N+NJET
 172   250   P(I,5)=SQRT(MAX(0.,P(I,4)**2-P(I,1)**2-P(I,2)**2-P(I,3)**2))
 173       ENDIF
 174
 175 C...Information about storage.
 176       MSTU(61)=N+1
 177       MSTU(62)=NP
 178       MSTU(63)=NC-N
 179       IF(MSTU(43).LE.1) MSTU(3)=NJET
 180       IF(MSTU(43).GE.2) N=N+NJET
 181
 182       RETURN
 183       END