TEvtGen/EvtGenModels/pykcut.F

   1
   2 C*********************************************************************
   3
   4 C...PYKCUT
   5 C...Dummy routine, which the user can replace in order to make cuts on
   6 C...the kinematics on the parton level before the matrix elements are
   7 C...evaluated and the event is generated. The cross-section estimates
   8 C...will automatically take these cuts into account, so the given
   9 C...values are for the allowed phase space region only. MCUT=0 means
  10 C...that the event has passed the cuts, MCUT=1 that it has failed.
  11
  12       SUBROUTINE PYKCUT(MCUT)
  13
  14 C...Double precision and integer declarations.
  15       IMPLICIT DOUBLE PRECISION(A-H, O-Z)
  16       IMPLICIT INTEGER(I-N)
  17       INTEGER PYK,PYCHGE,PYCOMP
  18 C...Commonblocks.
  19       COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
  20       COMMON/PYINT1/MINT(400),VINT(400)
  21       COMMON/PYINT2/ISET(500),KFPR(500,2),COEF(500,20),ICOL(40,4,2)
  22       COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
  23       SAVE /PYDAT1/,/PYINT1/,/PYINT2/,/PYPARS/
  24
  25 C...Set default value (accepting event) for MCUT.
  26       MCUT=0
  27
  28 C...Read out subprocess number.
  29       ISUB=MINT(1)
  30       ISTSB=ISET(ISUB)
  31
  32 C...Read out tau, y*, cos(theta), tau' (where defined, else =0).
  33       TAU=VINT(21)
  34       YST=VINT(22)
  35       CTH=0D0
  36       IF(ISTSB.EQ.2.OR.ISTSB.EQ.4) CTH=VINT(23)
  37       TAUP=0D0
  38       IF(ISTSB.GE.3.AND.ISTSB.LE.5) TAUP=VINT(26)
  39
  40 C...Calculate x_1, x_2, x_F.
  41       IF(ISTSB.LE.2.OR.ISTSB.GE.5) THEN
  42         X1=SQRT(TAU)*EXP(YST)
  43         X2=SQRT(TAU)*EXP(-YST)
  44       ELSE
  45         X1=SQRT(TAUP)*EXP(YST)
  46         X2=SQRT(TAUP)*EXP(-YST)
  47       ENDIF
  48       XF=X1-X2
  49
  50 C...Calculate shat, that, uhat, p_T^2.
  51       SHAT=TAU*VINT(2)
  52       SQM3=VINT(63)
  53       SQM4=VINT(64)
  54       RM3=SQM3/SHAT
  55       RM4=SQM4/SHAT
  56       BE34=SQRT(MAX(0D0,(1D0-RM3-RM4)**2-4D0*RM3*RM4))
  57       RPTS=4D0*VINT(71)**2/SHAT
  58       BE34L=SQRT(MAX(0D0,(1D0-RM3-RM4)**2-4D0*RM3*RM4-RPTS))
  59       RM34=2D0*RM3*RM4
  60       RSQM=1D0+RM34
  61       RTHM=(4D0*RM3*RM4+RPTS)/(1D0-RM3-RM4+BE34L)
  62       THAT=-0.5D0*SHAT*MAX(RTHM,1D0-RM3-RM4-BE34*CTH)
  63       UHAT=-0.5D0*SHAT*MAX(RTHM,1D0-RM3-RM4+BE34*CTH)
  64       PT2=MAX(VINT(71)**2,0.25D0*SHAT*BE34**2*(1D0-CTH**2))
  65
  66 C...Decisions by user to be put here.
  67 C      WRITE (*,*) "PYKCUT : ",SQRT(SHAT)
  68       IF(SHAT.LE.PARP(2)*PARP(2)) MCUT=1
  69 C...Stop program if this routine is ever called.
  70 C...You should not copy these lines to your own routine.
  71 C      WRITE(MSTU(11),5000)
  72 C      IF(PYR(0).LT.10D0) STOP
  73
  74 C...Format for error printout.
  75  5000 FORMAT(1X,'Error: you did not link your PYKCUT routine ',
  76      &'correctly.'/1X,'Dummy routine in PYTHIA file called instead.'/
  77      &1X,'Execution stopped!')
  78
  79       RETURN
  80       END