* * $Id: pregen.F,v 1.1.1.1 1996/01/11 14:14:41 mclareni Exp $ * * $Log: pregen.F,v $ * Revision 1.1.1.1 1996/01/11 14:14:41 mclareni * Cojets * * #include "cojets/pilot.h" SUBROUTINE PREGEN C ***************** C-- HANDLES PREPARATION OF TABLES FOR HARD PT GENERATION #if defined(CERNLIB_SINGLE) IMPLICIT REAL (A-H,O-Z) #endif #if defined(CERNLIB_DOUBLE) IMPLICIT DOUBLE PRECISION (A-H,O-Z) #endif #include "cojets/alqgen.inc" #include "cojets/bkwg.inc" #include "cojets/data1.inc" #include "cojets/decpar.inc" #include "cojets/entrev.inc" #include "cojets/event.inc" #include "cojets/flapre.inc" #include "cojets/forgen.inc" #include "cojets/fratab.inc" #include "cojets/inmat.inc" #include "cojets/intype.inc" #include "cojets/itapes.inc" #include "cojets/keybre.inc" #include "cojets/kfact.inc" #include "cojets/nevol.inc" #include "cojets/nflav.inc" #include "cojets/nkinpt.inc" #include "cojets/nleave.inc" #include "cojets/parq.inc" #include "cojets/qcds.inc" #include "cojets/stable.inc" #include "cojets/thrfla.inc" #include "cojets/tleave.inc" #include "cojets/tpretb.inc" #include "cojets/transl.inc" #include "cojets/ttends.inc" #include "cojets/tweigh.inc" DIMENSION WGIN(2) DATA ICALL/0/ IF(ICALL.GT.0) RETURN ICALL=1 C-- INITIALIZATION FLAPRE=1. CALL SBLOCK C-- SET DEFAULT SIGINL - FROM LOG(S)**2 FIT TO PBAR-P DATA AT ISR C-- (NUCL.PHYS. B262 (1985) 689) AND 900 GEV (UA5, Z.PHYSIK C 32 C-- (1986) 153) IF(SIGINL.EQ.0) SIGINL=36.66 *+(50.3-36.66)*(LOG(ECM)**2-LOG(62.3)**2) */(LOG(900.)**2-LOG(62.3)**2) C S=ECM**2 NFLAVS=NFLAV LGLUS=7 IVAL=100 ICHA(1,1)=IVAL ICHA(2,1)=IVAL ICHA(1,2)=LGLUS ICHA(2,2)=IVAL ICHA(1,3)=LGLUS ICHA(2,3)=LGLUS IWGV(1)=1 IWGV(2)=4 IWGV(3)=5 NEVENT=0 JBOOK=1 IF(PTMAX.EQ.0.) PTMAX=ECM*.45 IF(PTMAX.GT.ECM*.45) PTMAX=ECM*.45 ALLAM2=LOG(ALAMB**2) NALQ=127 IF(KPRHEV.EQ.1) PTMIN=MAX(REAL(QZEV*1.01),REAL(QMAS(KFRFLA))) PTI=MAX(PTMIN,PTMGE) ALQZM=LOG((QZEV/ALAMB)**2) ALQI=(LOG(PTI**2)-ALLAM2)/ALQZM IF(ALQI.LT.1.) GO TO 500 ALQF=(LOG(PTMAX**2)-ALLAM2)/ALQZM YF=LOG(ALQF) DALQ=(ALQF-ALQI)/FLOAT(NALQ) Q2VECT(1)=PTI**2 DO 5 I=1,NALQ I1=I+1 5 Q2VECT(I1)=EXP((ALQI+DALQ*FLOAT(I))*ALQZM+ALLAM2) C-- FORCED PT GENERATION -- INITIALIZATION IF(PTMGE.GT.PTMAX) GO TO 700 IPTWGT=0 8 CONTINUE C ZMCUT=1. ZLW=ZMCUT/ECM ZUP=1.-ZLW DO 3 L=1,6 AMHEV(L)=QMAS(L) AM2HEV(L)=(QMAS(L))**2 3 THRFLA(L)=4.*AM2HEV(L) XMIN=(3.*PTMIN**2-5.*PT2INT)/S IF(KPRHEV.EQ.1) XMIN=(THRFLA(KFRFLA)-5.*PT2INT)/S IF(KPRHEV.EQ.2) XMIN=(AM2HEV(KFRFLA)-5.*PT2INT)/S IF(PTMGE.GT.PTMIN) XMIN=(3.*PTMGE**2-5.*PT2INT)/S XMINJ=XMIN ALQFJ=ALQF YFJ=YF CALL INQCDS C ALAMBT=ALAMB NFLATR=NFLAV QSQMXT=ECM**2 CALL INIQCD C C-- MAIN LOOP C --------- DO 6 J1=1,3 DO 6 J2=1,5 DO 6 J3=1,NALQ 6 WGALQ(J3,J2,J1)=0. DO 7 N=1,NALQ 7 WGALQT(N)=0. C FRAT=0.1 IF(KPRHEV.GT.0) FRAT=0.3 FRATB=FRATAB FRAT=MAX(FRATB,FRAT) CALL TIMELF(TEXEC) TCASC=(1.-FRAT)*(TEXEC-TLEAVE)+TLEAVE IF(NSIGMA.GT.0) TCASC=TLEAVE NSIGMA=MAX(NSIGMA,NLEAVE*100) NCAS=0 DO 100 NC=1,NCASC CALL TIMELF(TLEFT) IF(TLEFT.LT.TCASC) GO TO 101 IF(NKINPT.GT.NSIGMA) GO TO 101 NCAS=NCAS+1 DO 100 ICH=1,3 C-- ICH=1 Q-Q OR Q-QB, 2 Q-G, 3 G-G C-- IWG=1 Q-Q OR Q-QB DIFFERENT FLAVORS, 2 Q-Q SAME FLAVORS, 3 Q-QB C-- SAME FLAVORS, 4 Q-G, 5 G-G CALL INPART(ICHA(1,ICH),IFLAIN(1),XIN(1),PXIN(1),PYIN(1),WGIN(1)) IF(IPBAR.EQ.1.AND.ICHA(1,ICH).NE.LGLUS) IFLAIN(1)=-IFLAIN(1) CALL INPART(ICHA(2,ICH),IFLAIN(2),XIN(2),PXIN(2),PYIN(2),WGIN(2)) IF(XIN(1)*XIN(2).LE.XMIN) GO TO 100 IWG=IWGV(ICH) IF(ICH.GT.1) GO TO 15 IF(IFLAIN(1).EQ.IFLAIN(2)) IWG=IWG+1 IF(IFLAIN(1).EQ.(-IFLAIN(2))) IWG=IWG+2 15 CONTINUE WGTCAS=WGIN(1)*WGIN(2) IF(ICH.EQ.2) WGTCAS=WGTCAS*2. CALL MATCH(IWG,WGTCAS) 100 CONTINUE 101 NCASC=NCAS IF(NCASC.EQ.0) GO TO 800 CALL TIMELF(TPRETB) TPRETB=TEXEC-TPRETB C C-- DETERMINE PROBABILITY OF MINBIAS EVENTS WGTOR=0. WGTHD1=0. DO 102 IW=1,5 DO 102 N=1,NALQ WGTHD1=WGTHD1+WGALQ(N,IW,1) IF(WGALQ(N,IW,3).LT.1.E-30) GO TO 102 WGTOR=WGTOR+WGALQ(N,IW,1)/WGALQ(N,IW,2)*WGALQ(N,IW,3) 102 CONTINUE FRHARD=WGTHD1/FLOAT(NCASC)*(.19733)**2*10./SIGINL*FACTK IF(FRHARD.GT..3) GO TO 900 PROPT=WGTOR PROPT=PROPT/FLOAT(NCASC) PROPT=PROPT*(.19733)**2*10./SIGINL*FACTK PROPTW=PROPT/(1.-(FRHARD-PROPT)) IF(PROPTW.GE.1.) GO TO 600 IF(PTMGE.LE.PTMIN) GO TO 105 C-- FORCED GENERATION -- HARD X-SECTION OVER FORCED PT RANGE SIGHRD=0. DO 104 IW=1,5 DO 104 N=1,NALQ 104 SIGHRD=SIGHRD+WGALQ(N,IW,1) SIGHRD=SIGHRD/FLOAT(NCASC)*(.19733)**2*10.*FACTK 105 CONTINUE C-- MODIFY PT TABLE WITH (OPTIONAL) PTWGT(PT) IF(IPTWGT.EQ.0) FRAGEH=PROPTW IF(IPTWGT.EQ.0) GO TO 110 WGTPT=0. DO 103 N=1,NALQ WPT=PTWGT(Q2VECT(N)) DO 103 IW=1,5 WGALQ(N,IW,1)=WGALQ(N,IW,1)*WPT WGALQ(N,IW,2)=WGALQ(N,IW,2)*WPT WGALQ(N,IW,3)=WGALQ(N,IW,3)*WPT 103 WGTPT=WGTPT+WGALQ(N,IW,1) PTWGTM=WGTPT/WGTHD1 110 CONTINUE TWEIGH=0. TENTRS=0. TWEIMN=0. C-- COMPLETE TABLES DO 120 N=1,NALQ WGALQT(N)=0. DO 120 IW=1,5 IF(WGALQ(N,IW,3).LT.1.E-30) GO TO 120 WGALQ3=WGALQ(N,IW,1) WGALQ(N,IW,2)=WGALQ(N,IW,1)/WGALQ(N,IW,2) WGALQ(N,IW,1)=WGALQ(N,IW,2)*WGALQ(N,IW,3) WGALQ(N,IW,3)=WGALQ3 WGALQT(N)=WGALQT(N)+WGALQ(N,IW,1) 120 CONTINUE DO 130 N=1,NALQ IF(WGALQT(N).LT.1.E-30) GO TO 130 WGALQ(N,1,1)=WGALQ(N,1,1)/WGALQT(N) DO 131 IW=2,5 WGALQ(N,IW,1)=WGALQ(N,IW,1)/WGALQT(N)+WGALQ(N,IW-1,1) 131 CONTINUE 130 CONTINUE DO 180 N=2,NALQ 180 WGALQT(N)=WGALQT(N)+WGALQT(N-1) IF(WGALQT(NALQ).EQ.0.) GO TO 920 DO 190 N=1,NALQ 190 WGALQT(N)=WGALQT(N)/WGALQT(NALQ) WGALQ3=0. DO 200 IW=1,5 DO 200 N=1,NALQ WGALQ3=WGALQ3+WGALQ(N,IW,3) 200 CONTINUE DO 201 IW=1,5 DO 201 N=1,NALQ 201 WGALQ(N,IW,3)=WGALQ(N,IW,3)/WGALQ3 NEVOL=0 RETURN C C-- ABNORMAL EXIT 500 IF(INTYPE.EQ.1) GO TO 510 WRITE(ITLIS,501) PTMIN,QZEV 501 FORMAT(1H1,7HPTMIN =,E11.4,31H GEV IS SMALLER THAN QZEV =, 1E11.4,49H GEV (INITIAL SCALE FOR PARTON DENSITY EVOLUTION) 2//1X,27HSET PTMIN GREATER THAN QZEV 3//' OR DECREASE QZEV (AND THEN MODIFY PARTON DENSITIES TO BE CON', 48HSISTENT) ) STOP 510 WRITE(ITLIS,511) PTMIN,QZEV 511 FORMAT(1H1,7HPTMIN =,E11.4,31H GEV IS SMALLER THAN QZEV =, 1E11.4,49H GEV (INITIAL SCALE FOR PARTON DENSITY EVOLUTION) 2//1X,'SET PTCUTOFF GREATER THAN QZEV' 3//' OR DECREASE QZEV (AND THEN MODIFY PARTON DENSITIES TO BE CON', 48HSISTENT) ) STOP 600 WRITE(ITLIS,601) PROPTW 601 FORMAT(1H1,22HQCD X-SECTION TOO HIGH * //1X,40HPROBABILITY FOR HARD EVENTS EXCEEDS 1 (= ,E10.3,1H) * //1X,47HINCREASE PTMIN (AND ALSO CHECK VALUE OF SIGINL) ) STOP 700 IF(INTYPE.EQ.1) GO TO 710 WRITE(ITLIS,701) PTMGE,PTMAX 701 FORMAT(1H1,14HINPUT PTMGE (= ,E10.3,26H) LARGER THAN MAX ALLOWED , *7HVALUE ( ,E10.3,1H) * //1X,14HDECREASE PTMGE ) STOP 710 WRITE(ITLIS,711) PTMGE,PTMAX 711 FORMAT(1H1,'INPUT PTFFFM (= ',E10.3,') LARGER THAN PTFFFX (' *,E10.3,1H) *//1X,'DECREASE PTFFFM') STOP 800 WRITE(ITLIS,801) 801 FORMAT(1H1,30HTOO LITTLE TIME FOR EXECUTION * //1X,36HINCREASE TIME PARAMETER ON JOB CARD ) STOP 900 IF(INTYPE.EQ.1) GO TO 910 WRITE(ITLIS,901) FRHARD,PTMIN 901 FORMAT(1H1,16HQCD X-SECTION = ,E12.4,7H*SIGINL * //1X,39HTOO HIGH (MAX LIMIT SET AT 0.3*SIGINL) * //1X,35HINCREASE PTMIN (AT PRESENT PTMIN = ,E12.3,1H) ) STOP 910 WRITE(ITLIS,911) FRHARD,PTMIN 911 FORMAT(1H1,16HQCD X-SECTION = ,E12.4,'*SIGTINEL' *//1X,'TOO HIGH (MAX LIMIT SET AT 0.3*SIGTINEL)' *//1X,'INCREASE PTCUTOFF (AT PRESENT PTCUTOFF = ',E12.3,1H) ) STOP 920 IF(INTYPE.EQ.1) GO TO 930 WRITE(ITLIS,921) 921 FORMAT(1H1,'NO ENTRIES IN PRETABULATION' *//1X,'TRY INCREASING NLEAVE, OR EXECUTION TIME AND/OR FRATAB') STOP 930 WRITE(ITLIS,931) 931 FORMAT(1H1,'NO ENTRIES IN PRETABULATION' */1X,'TRY INCREASING NUMBER OF REQUESTED EVENTS OR NSIGMA') STOP END