HepSim

Repository with Monte Carlo simulations for particle physics

jets_nlojetpp.py (raw text file)
# HepSim script using DataMelt http://jwork.org/dmelt
# Part of =HepMC= : https://atlaswww.hep.anl.gov/asc/hepsim/
# S.Chekanov (ANL)
from java.lang import *
from proto import FileNLO  # import FileNLO
from java.awt import Color,Font
from jhplot import  HPlot,P1D,HLabel,H1D
from jhplot.io import HBook
from jhplot.utils import FileList
from hephysics.particle import LParticle
import glob,os, sys


from javax.swing import JOptionPane
if JOptionPane.showConfirmDialog(None, "This macro works for local files. Quit?","Warning", JOptionPane.OK_CANCEL_OPTION) == JOptionPane.OK_OPTION: System.exit(0);


url=None
if len(sys.argv)>1: # if one argument is given, read all files from a directory 
   DIR=sys.argv[1]
   flist50=glob.glob(DIR+"*_pt50_*promc")
   flist100=glob.glob(DIR+"*_pt100_*promc")
   flist200=glob.glob(DIR+"*_pt200_*promc")
   flist400=glob.glob(DIR+"*_pt400_*promc")
   flist800=glob.glob(DIR+"*_pt800_*promc")
   flist1600=glob.glob(DIR+"*_pt1600_*promc")
   flist3200=glob.glob(DIR+"*_pt3200_*promc")
   print len(flist50),len(flist100),len(flist200),len(flist800),len(flist1600),len(flist3200)
   if (len(flist50)*len(flist100)*len(flist200)*len(flist800)*len(flist1600)*len(flist3200)==0):
               print "Error: some files are missing !"; sys.exit(0)
else:
   url="http://atlaswww1.hep.anl.gov/asc/RefHepSim/events"
   flist=[url+"/pp/8tev/jets_nlojetpp/nlojetpp_001.promc"]
   print "Reading 1 file from WWW=",flist[0]
   if len(flist)==0: print "Error: No input file!"; sys.exit(0)

minPT=50
maxETA=3.5
print "Inc. cross section. minPT=",minPT," max Eta=",maxETA

SYST=7+40 # 0, central, 7- scale uncertainties, 40-PDF (MSTW) 
bins=[50,100,200,400,800,1600,3200,6400]
# bins=[50,100,200,300,400,600,800,1200,1600,3200,6400]
# bins=[50,75,100,150,200,300,400,600,800,1200,1600,2400,3200,6400]
# bins=[50,75,100,150,200,300,400,600,800,1000,1200,1400,1600,2000,2400,3200,4000]

nev=0
ptbins=[[50,100], [100,200],[200,400], [400,800],[800,1600],[1600,3200],[3200,6400]]
flist= [ flist50, flist100, flist200,  flist400, flist800,   flist1600,  flist3200]

pthisto=[]
print "Process all pT ranges:",ptbins

for ipt in range(len(flist)):
  ptlist=flist[ipt]
  print "Dealing with="+str(ptbins[ipt])+" GeV"
  Nfiles=len(ptlist)
  totruns=0
  avcross=0
  nfiles=0
  histo=[]    # histograms with different PDF for direct
  for i in range(SYST):
      histo.append(  H1D("PT ="+str(i),bins) )

  for m in range(Nfiles):           # loop over all files in this list    
#    if (m>1): break; # debugging mode 
    bsize=os.path.getsize(ptlist[m])
    if (bsizeptbins[ipt][0] and pt

HEP.ANL.GOV