diff -Naur POWHEG-w/analize-W.f POWHEG-w-new/analize-W.f --- POWHEG-w/analize-W.f 2008-04-08 10:48:52.000000000 +0200 +++ POWHEG-w-new/analize-W.f 2009-03-24 15:02:38.000000000 +0100 @@ -240,8 +240,8 @@ real *8 cosths,deltaphi,deltaphi2 integer ifs,ihep integer maxtrack,maxjet - parameter (maxtrack=5000) - parameter (maxjet=100) + parameter (maxtrack=2048) + parameter (maxjet=2048) real *8 jetetcut,frenco,yktcut real *8 ktdeno,hacone,rparam diff -Naur POWHEG-w/bbarra.f POWHEG-w-new/bbarra.f --- POWHEG-w/bbarra.f 2008-05-19 18:05:19.000000000 +0200 +++ POWHEG-w-new/bbarra.f 2009-03-24 16:04:20.000000000 +0100 @@ -44,13 +44,14 @@ save firstcall,pborn,virt real *8 kt2,pt2plus,pt2minus,kt2l_damp external kt2,kt2l_damp - real *8 xr1,vr1,phir1,xr2,vr2,phir2 - common/cradvars2/xr1,vr1,phir1,xr2,vr2,phir2 + real *8 xr1,vr1,phir1,xr2,vr2,phir2,xjacr1,xjacr2 + common/cradvars2/xr1,vr1,phir1,xr2,vr2,phir2,xjacr1,xjacr2 cccccccccccccccccccccccccc c nlo histos variables real *8 bvweight,cweight,rweight,ctweight + save bvweight logical nloflag,outflag common/cnloflag/nloflag,outflag ccccccccccccccccccccccccc @@ -138,12 +139,19 @@ do j=-5,5 condition=abs(charge(i)+charge(j)-1d0).lt.2d-13 if(condition) then - bvweight=bvweight+(born(i,j)+virt(i,j)) + bvweight=bvweight+born(i,j) + if(order.ne.1) then + bvweight=bvweight+virt(i,j) + endif endif enddo enddo bvweight=bvweight*xjacb*hc2 - call outfun(w*bvweight,pborn) +c IMPORTANT: +c Do not call outfun here in order to fill NLO histograms +c since the weight w can vary during folded integration in MINT. +c The corresponding call will be made later, togheter with other +c Born-like configurations endif cccccccccccccccccc endif @@ -165,6 +173,17 @@ c The following lines are executed either if iflag=0 c or if iflag=1 + if(nloflag) then + if(.not.remflag) then + if (order.ne.3) then + call outfun(w*bvweight,pborn) +c The previous call is here only to comply with MINT with +c folded integration. See above for further details. + endif + endif + endif + + if(order.ne.1) then C REAL-LIKE EMISSION IN FIRST REGION OF SINGULARITY @@ -178,7 +197,8 @@ xr1=x vr1=v phir1=phi - + xjacr1=xjacrplus + if(check) then x_computed = 1-(dotp(preal(0,1),preal(0,5)) # +dotp(preal(0,2),preal(0,5)))/ @@ -222,6 +242,7 @@ cweight=cweight*xjacb*hc2 ctweight=ctweight*xjacb*hc2 + call outfun(w*cweight,pborn) call outfun(w*ctweight,pborn) @@ -268,6 +289,7 @@ xr2=x vr2=v phir2=phi + xjacr2=xjacrminus if(check) then x_computed = 1-(dotp(preal(0,1),preal(0,5))+ diff -Naur POWHEG-w/hevgen.f POWHEG-w-new/hevgen.f --- POWHEG-w/hevgen.f 2008-05-19 18:06:29.000000000 +0200 +++ POWHEG-w-new/hevgen.f 2009-03-24 17:18:03.000000000 +0100 @@ -789,8 +789,8 @@ common/cbornvars/ybar,Msq,costhlep,xpbar,xmbar logical remflag common/cremflag/remflag - real *8 xr1,vr1,phir1,xr2,vr2,phir2 - common/cradvars2/xr1,vr1,phir1,xr2,vr2,phir2 + real *8 xr1,vr1,phir1,xr2,vr2,phir2,xjacr1,xjacr2 + common/cradvars2/xr1,vr1,phir1,xr2,vr2,phir2,xjacr1,xjacr2 real *8 x,v,phi common/cradvars/x,v,phi real *8 x1,x2,pt2tmp1,pt2tmp2 @@ -842,9 +842,7 @@ call build_mom(ybar,Msq,costhlep,xr1,vr1,phir1,region,ptemp1 # ,x1,x2) -c At this point real matrix elements MUST be evaluated with -c alpha_s=alpha_s(pt^2) -C THIS CHANGE qfac TOO + x=xr1 v=vr1 phi=phir1 @@ -868,9 +866,6 @@ call build_mom(ybar,Msq,costhlep,xr2,vr2,phir2,region,ptemp2 # ,x1,x2) -c At this point real matrix elements MUST be evaluated with -c alpha_s=alpha_s(pt^2) -C THIS CHANGE qfac TOO x=xr2 v=vr2 phi=phir2 @@ -890,10 +885,13 @@ # +rminus2(jflborn(1),jflborn(2)) - rrr=random()*(valuep+valuem) - - if (rrr.lt.valuep) then + rrr=random()*(valuep*xjacr1+valuem*xjacr2) + if (rrr.lt.valuep*xjacr1) then +c The emission region must be chosen according to the full real +c contribution times the radiation phase space, which differs among +c the two regions. The following flavour choice (gq or qbarq) does +c not depend on it anymore. rrr2=random()*valuep if(rrr2.lt.rplus1(jflborn(1),jflborn(2))) then diff -Naur POWHEG-w/lhefwrite.f POWHEG-w-new/lhefwrite.f --- POWHEG-w/lhefwrite.f 2008-01-29 11:29:58.000000000 +0100 +++ POWHEG-w-new/lhefwrite.f 2009-03-24 15:03:12.000000000 +0100 @@ -17,7 +17,7 @@ call rm48ut(iran,n1ran,n2ran) write(nlf,*) 'Random number generator initialized with: ', # iran,' ',n1ran,' ',n2ran - write(nlf,'(a)') '!-->' + write(nlf,'(a)') '-->' write(nlf,'(a)') '' write(nlf,110) idbmup(1),idbmup(2),ebmup(1),ebmup(2), &pdfgup(1),pdfgup(2),pdfsup(1),pdfsup(2),idwtup,nprup diff -Naur POWHEG-w/main-HERWIG.f POWHEG-w-new/main-HERWIG.f --- POWHEG-w/main-HERWIG.f 2009-03-24 17:59:39.000000000 +0100 +++ POWHEG-w-new/main-HERWIG.f 2009-03-24 19:27:33.000000000 +0100 @@ -2,16 +2,22 @@ C---COMMON BLOCKS ARE INCLUDED AS FILE herwig6510.h INCLUDE 'herwig6510.h' include 'LesHouches.h' - real * 8 powheginput - EXTERNAL powheginput INTEGER N logical uevent parameter (uevent=.true.) + character * 20 pwgprefix + integer lprefix + common/cpwgprefix/pwgprefix,lprefix + real * 8 powheginput + external powheginput C---PROCESS; set to negative for user supplied me iproc=-1 ! Les Houches interface C---MAX NUMBER OF EVENTS IN THIS RUN; must be set before HWIGIN call maxev=powheginput('maxev') + open(unit=23,file=pwgprefix(1:lprefix)//'uboundfailures.dat', + # status='unknown') + write(23,*) 'pt, mu, value, ubound, upper' C---INITIALISE OTHER COMMON BLOCKS CALL HWIGIN C---USER CAN RESET PARAMETERS AT @@ -103,4 +109,12 @@ call aend close(99) call printstat +c save last random number + open(unit=24,file=pwgprefix(1:lprefix)//'rm48.out', + # status='unknown') + call rm48ut(iran,n1ran,n2ran) + write(24,*) 'Random number generator values to resume this run :', + & iran,' ',n1ran,' ',n2ran + close(24) + close(23) end diff -Naur POWHEG-w/main-lhef.f POWHEG-w-new/main-lhef.f --- POWHEG-w/main-lhef.f 2009-03-24 17:59:39.000000000 +0100 +++ POWHEG-w-new/main-lhef.f 2009-03-24 19:27:11.000000000 +0100 @@ -42,7 +42,7 @@ open(unit=24,file=pwgprefix(1:lprefix)//'rm48.out', # status='unknown') call rm48ut(iran,n1ran,n2ran) - write(24,*) 'Random number generator values to resum this run :', + write(24,*) 'Random number generator values to resume this run :', & iran,' ',n1ran,' ',n2ran close(24) close(17) diff -Naur POWHEG-w/main-PYTHIA.f POWHEG-w-new/main-PYTHIA.f --- POWHEG-w/main-PYTHIA.f 2009-03-24 17:59:39.000000000 +0100 +++ POWHEG-w-new/main-PYTHIA.f 2009-03-24 19:29:56.000000000 +0100 @@ -18,10 +18,12 @@ & VHEP(4,NMXHEP) integer iev,temp + external pydata + character * 20 pwgprefix + integer lprefix + common/cpwgprefix/pwgprefix,lprefix real * 8 powheginput external powheginput - external pydata - c PARAMETERS mstj(41)=3 !No photon radiations off leptons c mstp(61)=0 !No IS shower @@ -32,9 +34,12 @@ c mstp(131)=0 !No Pile Up c mstp(111)=0 !No hadronization + open(unit=23,file=pwgprefix(1:lprefix)//'uboundfailures.dat', + # status='unknown') + write(23,*) + # 'ifl,jfl, pt, mu, value, ubound, region, lum_fac, upper' - -c Set up PYTHIA to accept user processes +c Set up PYTHIA to accept user processes call PYINIT('USER','','',0d0) maxev=powheginput('maxev') call PYABEG @@ -69,6 +74,14 @@ write(*,*) 'At the end NEVHEP is ',nevhep call PYAEND call printstat +c save last random number + open(unit=24,file=pwgprefix(1:lprefix)//'rm48.out', + # status='unknown') + call rm48ut(iran,n1ran,n2ran) + write(24,*) 'Random number generator values to resume this run :', + & iran,' ',n1ran,' ',n2ran + close(24) + close(23) END subroutine UPINIT diff -Naur POWHEG-w/Makefile POWHEG-w-new/Makefile --- POWHEG-w/Makefile 2009-03-24 18:01:10.000000000 +0100 +++ POWHEG-w-new/Makefile 2009-03-24 17:40:03.000000000 +0100 @@ -11,8 +11,8 @@ #F77= gfortran -Wall -fno-automatic -ffpe-trap=invalid,zero,overflow # -fbounds-check # The following should no longer be needed in gfortran >= 4.3 #TRAPFPE= -F77= g77 -Wall -fno-automatic -#-ffortran-bounds-check +F77= g77 -Wall +#-fno-automatic -ffortran-bounds-check DEBUG= #-ggdb -pg TRAPFPE=trapfpe.o diff -Naur POWHEG-w/mint-integrator.f POWHEG-w-new/mint-integrator.f --- POWHEG-w/mint-integrator.f 2008-01-29 11:29:58.000000000 +0100 +++ POWHEG-w-new/mint-integrator.f 2009-03-24 17:06:40.000000000 +0100 @@ -184,7 +184,7 @@ elseif(err.eq.0) then err=etot endif - ans=(ans/err+vtot/etot)/(1/err+1/etot) + ans=(ans/err**2+vtot/etot**2)/(1/err**2+1/etot**2) err=1/sqrt(1/err**2+1/etot**2) endif goto 10 diff -Naur POWHEG-w/POWHEG-w.f POWHEG-w-new/POWHEG-w.f --- POWHEG-w/POWHEG-w.f 2009-03-24 17:59:39.000000000 +0100 +++ POWHEG-w-new/POWHEG-w.f 2009-03-24 15:37:18.000000000 +0100 @@ -147,7 +147,8 @@ real *8 alfas external alfas common/chardrad/preal,pt2val - + real *8 tiny + parameter(tiny=1d-5) call choosegen @@ -157,7 +158,7 @@ scalup=sqrt(max(pt2val,ptminsq)) c Set up the Les Houches common block - if(pt2val.lt.ptminsq/2) then + if(pt2val.lt.tiny) then c meaning: pt2val=0!!! NUP=3 call set_ren_fac_scales ! needed only to write correct aqcdup