TEvtGen/EvtGenModels/jetset1.F

   1 C--------------------------------------------------------------------------
   2 C
   3 C Environment:
   4 C      This software is part of the EvtGen package developed jointly
   5 C      for the BaBar and CLEO collaborations.  If you use all or part
   6 C      of it, please give an appropriate acknowledgement.
   7 C
   8 C Copyright Information: See EvtGen/COPYRIGHT
   9 C      Copyright (C) 1998      Caltech, UCSB
  10 C
  11 C Module: jetset1.F
  12 C
  13 C Description:
  14 C
  15 C Modification history:
  16 C
  17 C    DJL/RYD     August 11, 1998         Module created
  18 C
  19 C------------------------------------------------------------------------
  20       subroutine jetset1(ip,m,ndaug,kf,km,px,py,pz,e)
  21
  22 C
  23 C interface to JETSET 7.4 to have one particle decayed
  24 C including possibly fragmentation, if the decay products include
  25 C partons.
  26 C
  27
  28       implicit none
  29
  30       INTEGER MSTU,MSTJ
  31       REAL PARU,PARJ
  32       COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
  33       SAVE /LUDAT1/
  34
  35       common/lujets/n,k(4000,5),p(4000,5),v(4000,5)
  36
  37       integer n,k
  38       real p,v
  39
  40       integer ndmax
  41       parameter (ndmax=100)
  42
  43       integer lucomp
  44       external lucomp
  45
  46       integer ip,kf(ndmax),i,ndaug,km(ndmax)
  47       integer kp,kid,ipart1,ipart
  48       real*8  m,px(ndmax),py(ndmax),pz(ndmax),e(ndmax)
  49
  50       integer ip1
  51       real qmax
  52
  53 c      used to use lu1ent both since it does not set the mass
  54 c      of the daughter we have to manipulate the common blocks our
  55 c      self. ryd April 25-1999
  56
  57 c      call lu1ent(1,ip,0.0,0.0,0.0)
  58
  59       K(1,1)=1
  60       K(1,2)=ip
  61       P(1,5)=m
  62       P(1,4)=m
  63       P(1,1)=0.0
  64       P(1,2)=0.0
  65       P(1,3)=0.0
  66       n=1
  67
  68
  69 c     now we can decay this particle
  70       call ludecy(1)
  71
  72 C code copied from LUEXEC to avoid error with shower
  73 C switched on
  74 C...Decay products may develop a shower.
  75         IF(MSTJ(92).GT.0) THEN
  76           IP1=MSTJ(92)
  77           QMAX=SQRT(MAX(0.,(P(IP1,4)+P(IP1+1,4))**2-(P(IP1,1)+P(IP1+1,
  78      &    1))**2-(P(IP1,2)+P(IP1+1,2))**2-(P(IP1,3)+P(IP1+1,3))**2))
  79           CALL LUSHOW(IP1,IP1+1,QMAX)
  80           CALL LUPREP(IP1)
  81           MSTJ(92)=0
  82         ELSEIF(MSTJ(92).LT.0) THEN
  83           IP1=-MSTJ(92)
  84           CALL LUSHOW(IP1,-3,P(IP1,5))
  85           CALL LUPREP(IP1)
  86           MSTJ(92)=0
  87         ENDIF
  88
  89 c
  90 c for debugging:
  91 c      call lulist(1)
  92 c
  93       mstj(21)=0
  94       call luexec
  95       mstj(21)=2
  96
  97 c find partons, delete secondary partons, set mother pointers
  98
  99       ndaug = 0
 100       ipart1 = 1
 101       ipart = 1
 102
 103       do 10 i=2,n
 104         kp = k(i,3)
 105         kid = k(i,2)
 106         if (abs(kid) .ge. 1 .and. abs(kid) .le. 8
 107      1  .or. kid .eq. 21
 108      2  .or. kid .ge. 91 .and. kid .le. 94)  then
 109            if (ipart1 .eq. 1) ipart1 = i
 110            ipart = i
 111            if (kp .ne. 1)  goto 10
 112            kp = 0
 113         else
 114            if (kp .gt. ipart)  then
 115               goto 10
 116            elseif (kp .ge. ipart1)  then
 117               kp = ipart1-1
 118            else
 119               kp = 0
 120               endif
 121            endif
 122         ndaug = ndaug + 1
 123         km(ndaug)=kp
 124         kf(ndaug)=kid
 125         px(ndaug)=p(i,1)
 126         py(ndaug)=p(i,2)
 127         pz(ndaug)=p(i,3)
 128         e(ndaug)=p(i,4)
 129 c
 130 c for debugging:
 131 c        print '( 2I5,I12,4F12.4 )',ndaug,km(ndaug),kf(ndaug),
 132 c     1    px(ndaug),py(ndaug),pz(ndaug),e(ndaug)
 133 c
 134    10 continue
 135
 136       end
 137
 138