C----------------------------------------------------------------------
C
C Code to calculate the PDFs
C
C PDF actaully calculates the PDFs at scale EMSCAB and x_1=XXB(1)
C x_2=XXB(2)
C is returns DISFB(13,2) where the second index is the beam particle
C and the first is the type of parton
C 1, 2, 3, 4, 5, 6 = d,u,s,b,c,b,t
C 7, 8, 9,10,11,12 = dbar,ubar,sbar,cbar,bbar,tbar
C 13=gluon
C
C there are also dummy subroutines to use if PDFLIB is not linked.
C The code may be specified with the following options
C
C -DPYTHIA use PYTHIA default PDFs
C -DHERWIG use HERWIG default PDFs
C -DPDFLIB use PDFLIB
C
C in each case the correct code must be linked.
C
C----------------------------------------------------------------------
CDECK ID>, PDF.
*-- Author : Peter Richardson
C----------------------------------------------------------------------
SUBROUTINE PDF(EMSCAB,XXB,DISFB)
C----------------------------------------------------------------------
C Subroutine to calculate the pdfs
C generator dependent
C----------------------------------------------------------------------
#if defined(HERWIG)
INCLUDE 'HERWIG65.INC'
DOUBLE PRECISION EMSCAB,XXB(2),DISFB(13,2),I,J,
& BHRUNI,UPV,DNV,USEA,DSEA,STR,CHM,BTM,TOP,GLU
INTEGER IDB(2)
CHARACTER *8 DUMMY
DOUBLE PRECISION VALUE(20)
CHARACTER*20 PARM(20)
LOGICAL FIRST
EXTERNAL BHRUNI
DATA VALUE/20*0D0/,PARM/20*' '/
DATA FIRST/.TRUE./
C--Les Houches run common block
INTEGER MAXPUP
PARAMETER(MAXPUP=100)
INTEGER IDBMUP,PDFGUP,PDFSUP,IDWTUP,NPRUP,LPRUP
DOUBLE PRECISION EBMUP,XSECUP,XERRUP,XMAXUP
COMMON /HEPRUP/ IDBMUP(2),EBMUP(2),PDFGUP(2),PDFSUP(2),
& IDWTUP,NPRUP,XSECUP(MAXPUP),XERRUP(MAXPUP),
& XMAXUP(MAXPUP),LPRUP(MAXPUP)
C--use generator internal pdfs
IF(PDFGUP(1).LT.0) THEN
C---SET UP INITIAL STATE
CALL HWUIDT(1,IDBMUP(1),IDB(1),DUMMY)
PART1=DUMMY
CALL HWUIDT(1,IDBMUP(2),IDB(2),DUMMY)
PART2=DUMMY
IPART1 = IDB(1)
IPART2 = IDB(2)
NHEP=1
ISTHEP(NHEP)=101
PHEP(1,NHEP)=0.
PHEP(2,NHEP)=0.
PHEP(3,NHEP)=EBMUP(1)
PHEP(4,NHEP)=EBMUP(1)
PHEP(5,NHEP)=RMASS(IDB(1))
JMOHEP(1,NHEP)=0
JMOHEP(2,NHEP)=0
JDAHEP(1,NHEP)=0
JDAHEP(2,NHEP)=0
IDHW(NHEP)=IPART1
IDHEP(NHEP)=IDPDG(IDB(1))
NHEP=NHEP+1
ISTHEP(NHEP)=102
PHEP(1,NHEP)=0.
PHEP(2,NHEP)=0.
PHEP(3,NHEP)=-EBMUP(2)
PHEP(4,NHEP)=EBMUP(2)
PHEP(5,NHEP)=RMASS(IDB(2))
JMOHEP(1,NHEP)=0
JMOHEP(2,NHEP)=0
JDAHEP(1,NHEP)=0
JDAHEP(2,NHEP)=0
IDHW(NHEP)=IPART2
IDHEP(NHEP)=IDPDG(IPART2)
C---NEXT ENTRY IS OVERALL CM FRAME
NHEP=NHEP+1
IDHW(NHEP)=14
IDHEP(NHEP)=0
ISTHEP(NHEP)=103
JMOHEP(1,NHEP)=NHEP-2
JMOHEP(2,NHEP)=NHEP-1
JDAHEP(1,NHEP)=0
JDAHEP(2,NHEP)=0
CALL BHVSUM(4,PHEP(1,NHEP-1),PHEP(1,NHEP-2),PHEP(1,NHEP))
CALL BHUMAS(PHEP(1,NHEP))
XLMIN = XLMINB
EMSCA = EMSCAB
DO I=1,2
XX(I)=XXB(I)
ENDDO
CALL HWSGEN(.TRUE.)
DO I=1,2
DO J=1,13
DISFB(J,I) = DISF(J,I)
ENDDO
ENDDO
ELSE
#endif
#if defined(PYTHIA)
IMPLICIT NONE
C--Les Houches run common block
INTEGER MAXPUP
PARAMETER(MAXPUP=100)
INTEGER IDBMUP,PDFGUP,PDFSUP,IDWTUP,NPRUP,LPRUP
DOUBLE PRECISION EBMUP,XSECUP,XERRUP,XMAXUP
COMMON /HEPRUP/ IDBMUP(2),EBMUP(2),PDFGUP(2),PDFSUP(2),
& IDWTUP,NPRUP,XSECUP(MAXPUP),XERRUP(MAXPUP),
& XMAXUP(MAXPUP),LPRUP(MAXPUP)
DOUBLE PRECISION EMSCAB,XXB(2),DISFB(13,2),I,J,XXMINB,
& XLMINB,BHRUNI,XPQ(-25:25),Q2,UPV,DNV,USEA,DSEA,STR,CHM,BTM,
& TOP,GLU
EXTERNAL BHRUNI
DOUBLE PRECISION VALUE(20)
CHARACTER*20 PARM(20)
LOGICAL FIRST
DATA VALUE/20*0D0/,PARM/20*' '/
DATA FIRST/.TRUE./
IF(PDFGUP(1).LT.0) THEN
XXB(1)=EXP(BHRUNI(0,0.0D0,XLMINB))
XXB(2)=XXMINB/XXB(1)
Q2 = EMSCAB**2
CALL PYPDFU(IDBMUP(1),XXB(1),Q2,XPQ)
DO I=1,6
DISFB(I,1) = XPQ(I)
DISFB(I+6,1) = XPQ(I-6)
ENDDO
DISFB(13,1) = XPQ(21)
CALL PYPDFU(IDBMUP(2),XXB(2),Q2,XPQ)
DO I=1,6
DISFB(I,2) = XPQ(I)
DISFB(I+6,2) = XPQ(I-6)
ENDDO
DISFB(13,2) = XPQ(21)
ELSE
#else
IMPLICIT NONE
C--Les Houches run common block
INTEGER MAXPUP
PARAMETER(MAXPUP=100)
INTEGER IDBMUP,PDFGUP,PDFSUP,IDWTUP,NPRUP,LPRUP
DOUBLE PRECISION EBMUP,XSECUP,XERRUP,XMAXUP
COMMON /HEPRUP/ IDBMUP(2),EBMUP(2),PDFGUP(2),PDFSUP(2),
& IDWTUP,NPRUP,XSECUP(MAXPUP),XERRUP(MAXPUP),
& XMAXUP(MAXPUP),LPRUP(MAXPUP)
DOUBLE PRECISION EMSCAB,XXB(2),DISFB(13,2),UPV,DNV,USEA,DSEA,
& STR,CHM,BTM,TOP,GLU,VALUE(20)
INTEGER I
CHARACTER*20 PARM(20)
LOGICAL FIRST
DATA VALUE/20*0D0/,PARM/20*' '/
DATA FIRST/.TRUE./
IF(PDFGUP(1).LT.0) THEN
print *,'testing'
stop
ELSE
#endif
C--initialisation
IF(FIRST) THEN
IF(PDFGUP(1).NE.PDFGUP(2).OR.PDFSUP(1).NE.PDFSUP(2)) THEN
print *,'MUST HAVE SAME PDF FOR BOTH BEAMS'
STOP
ENDIF
PARM(1)='NPTYPE'
VALUE(1)=1
PARM(2)='NGROUP'
VALUE(2)=PDFGUP(1)
PARM(3)='NSET'
VALUE(3)=PDFSUP(1)
CALL PDFSET(PARM,VALUE)
FIRST = .FALSE.
ENDIF
DO I=1,2
CALL STRUCTM(XXB(I),EMSCAB,
& UPV,DNV,USEA,DSEA,STR,CHM,BTM,TOP,GLU)
DISFB( 3,I) = STR
DISFB( 4,I) = CHM
DISFB( 5,I) = BTM
DISFB( 6,I) = TOP
DISFB( 9,I) = STR
DISFB(10,I) = CHM
DISFB(11,I) = BTM
DISFB(12,I) = TOP
DISFB(13,I) = GLU
IF(IDBMUP(I).GT.0) THEN
DISFB(1,I) = DNV+DSEA
DISFB(2,I) = UPV+USEA
DISFB(7,I) = DSEA
DISFB(8,I) = USEA
ELSE
DISFB(1,I) = DSEA
DISFB(2,I) = USEA
DISFB(7,I) = DNV+DSEA
DISFB(8,I) = UPV+USEA
ENDIF
ENDDO
ENDIF
END
#if !defined(PDFLIB)
CDECK ID>, PDFSET.
*CMZ :- -26/04/91 11.11.54 by Bryan Webber
*-- Author : Bryan Webber
C----------------------------------------------------------------------
SUBROUTINE PDFSET(PARM,VAL)
C----------------------------------------------------------------------
C DUMMY SUBROUTINE: DELETE AND SET MODPDF(I)
C IN MAIN PROGRAM IF YOU USE PDFLIB CERN-LIBRARY
C PACKAGE FOR NUCLEON STRUCTURE FUNCTIONS
C----------------------------------------------------------------------
DOUBLE PRECISION VAL(20)
CHARACTER*20 PARM(20)
WRITE (6,10)
10 FORMAT(/10X,'PDFSET CALLED BUT NOT LINKED')
STOP
END
CDECK ID>, STRUCTM.
*CMZ :- -26/04/91 11.11.54 by Bryan Webber
*-- Author : Bryan Webber
C-----------------------------------------------------------------------
SUBROUTINE STRUCTM(X,QSCA,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU)
C-----------------------------------------------------------------------
C DUMMY SUBROUTINE: DELETE IF YOU USE PDFLIB CERN-LIBRARY
C PACKAGE FOR NUCLEON STRUCTURE FUNCTIONS
C-----------------------------------------------------------------------
DOUBLE PRECISION X,QSCA,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU
WRITE (6,10)
10 FORMAT(/10X,'STRUCTM CALLED BUT NOT LINKED')
STOP
END
#endif