C----------------------------------------------------------------------- PROGRAM ZGEN C----------------------------------------------------------------------- C Example main program C----------------------------------------------------------------------- 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) C--common block to pass the parameters DOUBLE PRECISION SMIN,SMAX INTEGER IPR COMMON /LPARAM/SMIN,SMAX,IPR C--local variables INTEGER I,NPOINT DOUBLE PRECISION WGT,WGTSM,WGTSQ PARAMETER(NPOINT=100) C--beam particles and energies (only proton=2212 and antiproton=-2212) IDBMUP(1) = 2212 IDBMUP(2) = 2212 EBMUP(1) = 7000.0D0 EBMUP(2) = 7000.0D0 C--pdf's for the beams (use CTEQ5L) C--MUST BE THE SAME FOR BOTH BEAMS PDFGUP(1) = 4 PDFGUP(2) = 4 PDFSUP(1) = 46 PDFSUP(2) = 46 C--integration parameters SMIN = 10.0D0**2 SMAX = (EBMUP(1)+EBMUP(2))**2 print *,'limits on sqrt(sh)',sqrt(smin),sqrt(smax) C--which process the following are supported C 100 all fermion-antifermion production C 110 all lepton production C 120 all quark production C 130+I production of the fermion i=(1,6 =e,ve,mu,vmu,tau,vtau) C (7,12-d,u,s,c,b,t) IPR = 110 C--set up the parameters of the Little Higgs Model CALL LITTLE C--code to compute the average weight WGTSM = 0.0D0 WGTSQ = 0.0D0 DO I=1,NPOINT CALL ZWGHT(WGT) WGTSM = WGTSM+WGT WGTSQ = WGTSQ+WGT**2 ENDDO WGTSM = WGTSM/DBLE(NPOINT) WGTSQ = MAX(0.0D0,WGTSQ/DBLE(NPOINT)-WGTSM**2) WGTSQ = SQRT(WGTSQ/DBLE(NPOINT)) print *,'Cross section = ',wgtsm,'+/-',wgtsq,' pb' END C----------------------------------------------------------------------- SUBROUTINE ZWGHT(WGT) C----------------------------------------------------------------------- C Subroutine to return the weight for a given phase space point C----------------------------------------------------------------------- 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) C--common block to pass the parameters DOUBLE PRECISION SMIN,SMAX INTEGER IPR COMMON /LPARAM/SMIN,SMAX,IPR C--little Higgs common block C--common block with Liggle Higgs parameters DOUBLE PRECISION MW,MZ,MF,ALPHA,SW,MT,FH,THETA,THETAP,G,GP,E,PI, & CW,VEV,VEVP,ST,CT,STP,CTP,LAM,MZH,MWH,MAH,GL,GR,WIDTH COMMON /LHIGGS/MW,MZ,MF(12),ALPHA,SW,MT,FH,THETA,THETAP,G,GP,E,PI, & CW,VEV,VEVP,ST,CT,STP,CTP, & LAM(2),MZH,MWH,MAH,GL(14,6),GR(14,6),WIDTH(6) C--local variables INTEGER NCH PARAMETER(NCH=4) INTEGER ITYPE(NCH),I,IFERM(3) DOUBLE PRECISION MCH(NCH),MCH2(NCH),MGAM(NCH),MGAM2(NCH),ALP(NCH), & WI(NCH),LIM(2,NCH),RANGEN DOUBLE PRECISION WGT,CTHETA,PHI,PIFAC,RANUNI,GEV2PB,TEMP,SH,XX(2), & FJAC,DISF(13,2) EXTERNAL RANGEN,RANUNI LOGICAL FIRST DATA FIRST/.TRUE./ SAVE FIRST,PIFAC,ITYPE,WI,ALP,MCH,MCH2,MGAM,MGAM2,LIM, & GEV2PB,IFERM C--useful parameters and initilisation IF(FIRST) THEN C--constants PIFAC = ACOS(-1.0D0) FIRST = .FALSE. C--conversion to picobarn from GeV GEV2PB=389379.D3 C--which processes to generate (IFERM controls the loop over final state fermions) C--(e,ve,mu,vmu,tau,nvtau,d,u,s,c,b,t) IF(IPR.EQ.100) THEN IFERM(1) = 1 IFERM(2) = 12 ELSEIF(IPR.EQ.110) THEN IFERM(1) = 1 IFERM(2) = 6 ELSEIF(IPR.EQ.120) THEN IFERM(1) = 7 IFERM(2) = 12 ELSEIF(IPR.GE.131.AND.IPR.LE.142) THEN IFERM(1) = IPR-130 IFERM(2) = IPR-130 ELSE print *,'unknown process stopping' STOP ENDIF print *,'selected process code =',ipr print *,'limits on fermion loop',iferm C--initialisation of the multichannel integration C--first a channel to smooth the photon ITYPE(1) = 1 WI (1) = 0.4D0 ALP (1) = 2.0D0 MCH (1) = 0.0D0 MCH2 (1) = 0.0D0 MGAM (1) = 0.0D0 MGAM2(1) = 0.0D0 C--second a channel to smooth the Z ITYPE(2) = 2 WI (2) = 0.4D0 ALP (2) = 0.0D0 MCH (2) = MZ MCH2 (2) = MZ**2 MGAM (2) = MZ*WIDTH(2) MGAM2(2) = MGAM(2)**2 C--third a channel to smooth the heavy photon ITYPE(3) = 2 WI (3) = 0.1D0 ALP (3) = 0.0D0 MCH (3) = MAH MCH2 (3) = MAH**2 MGAM (3) = MAH*WIDTH(3) MGAM2(3) = MGAM(3)**2 C--fourth a channel to smooth the heavy Z boson ITYPE(4) = 2 WI (4) = 0.1D0 ALP (4) = 0.0D0 MCH (4) = MZH MCH2 (4) = MZH**2 MGAM (4) = MZH*WIDTH(3) MGAM2(4) = MGAM(4)**2 C--compute the limits on the variables CTHETA = 0.0D0 DO I=1,NCH CTHETA = CTHETA+WI(I) IF(ITYPE(I).EQ.1) THEN LIM(1,I) = SMAX**(1.0D0-ALP(I)) LIM(2,I) = SMIN**(1.0D0-ALP(I)) ELSEIF(ITYPE(I).EQ.2) THEN LIM(1,I) = ATAN((SMIN-MCH2(I))/MGAM(I)) LIM(2,I) = ATAN((SMAX-MCH2(I))/MGAM(I)) ENDIF ENDDO C--ensure weights add up to one CTHETA = 1.0D0/CTHETA DO I=1,NCH WI(I)=WI(I)*CTHETA ENDDO ENDIF C--calculation of the phase space weight and momenta C--first the values of theta and phi CTHETA = RANUNI(1,-1.0D0,1.0D0) PHI = RANUNI(2,0.0D0,2.0D0*PIFAC) C--weight for this piece of the process WGT = 4.0D0*PIFAC C--calculation of sh C--first pick the channel TEMP = RANGEN(3) DO I=1,NCH IF(WI(I).GT.TEMP) GOTO 10 TEMP = TEMP-WI(I) ENDDO C--generate the transformed variable 10 TEMP = RANUNI(4,LIM(1,I),LIM(2,I)) C--calculate sh from transformed varaible IF(ITYPE(I).EQ.1) THEN SH = TEMP**(1.0D0/(1.0D0-ALP(I))) ELSE SH = MGAM(I)*TAN(TEMP)+MCH2(I) ENDIF C--compute the Jacobian function FJAC = 0.0D0 DO I=1,NCH TEMP = LIM(2,I)-LIM(1,I) IF(ITYPE(I).EQ.1) THEN TEMP = 1.0D0/SH**ALP(I)/TEMP ELSE TEMP = MGAM(I)/TEMP/((SH-MCH2(I))**2+MGAM2(I)) ENDIF FJAC = FJAC+WI(I)*TEMP ENDDO C--include this in the phase-space weight WGT = WGT/FJAC C--generate the value of x_1 and x_2 TEMP = LOG(SH/SMAX) XX(1) = EXP(RANUNI(5,LOG(SH/SMAX),0.0D0)) XX(2) = SH/SMAX/XX(1) C--include this in the phase-space weight WGT = -WGT*TEMP C--final factors for the phase space piece of the weight WGT = WGT/32.0D0/PIFAC**2/SH**2.5D0 C--call the structure functions CALL PDF(SQRT(SH),XX,DISF) print *,'testing the weight and sqrt(sh)',wgt,sqrt(sh) END