In [1]:
import lhef
import lorentz
In [2]:
import matplotlib.pyplot as plt
In [3]:
import plotly.plotly as py
In [4]:
from bokeh.plotting import figure, show, output_notebook, output_file, reset_output, save#, push_notebook
output_notebook()
In [5]:
import numpy as np
import math
import utils as u
import pandas as pd
In [6]:
import importlib
In [4]:
importlib.reload(lorentz)
Out[4]:
In [5]:
importlib.reload(lhef)
Out[5]:
In [7]:
def arange(a,b,s):
return np.arange(a,b+s,s)
In [262]:
LHEfile=lhef.readLHE("pp2zz24l.lhe")
In [269]:
LHEfile=lhef.readLHE("pp2h24l.lhe")
A first example: collect as you go.¶
In [223]:
nprinted=0
debug=True
nPrint=11
costheta_values=[] # create a vector where to store the computed values of costheta
abs_costheta_values=[] # and one for the abs
theta_values=[] # and one for the theta angle
for e in LHEfile: # loop on the events
for p in e.particles: # loop on the particles of each event
if p.status == 1 and p.id == 11: # check it is a final state and is a e-
lv=p.fourvector() # make four vector
obs=lv.cosTheta() # obtain the cosTheta
costheta_values.append(obs) # append it to the vector of results
abs_costheta_values.append(math.fabs(obs)) # and the abs(cosTheta)
theta=lv.theta() # obtain the theta angle
theta_values.append(theta)
if nprinted <nPrint:
if debug: print(p.px, p.py, p.pz, p.e)
if debug: print(obs)
nprinted+=1
In [236]:
[1,2,3,4][::3]
Out[236]:
In [229]:
[1,2,3,4][1:]
Out[229]:
In [225]:
hist, edges = np.histogram(costheta_values,bins=arange(-1,1,0.2),normed=True)
p1 = figure(title="Electron cosθ Distribution in BG",#tools="save",
background_fill_color="#E8DDCB")
p1.quad(top=hist, bottom=0, left=edges[:-1], right=edges[1:],
fill_color="#036564", line_color="#033649")
p1.xaxis.axis_label = 'cosθ'
p1.yaxis.axis_label = 'dσ/dcosθ'
show(p1)
In [129]:
plt.hist(costheta_values,bins=arange(-1,1,0.1),density=True)
plt.title("Electron θ Distribution")
plt.xlabel("cosθ")
plt.ylabel("dσ/dcosθ")
plt.show()
A second example: collect pairs as you go.¶
In [247]:
nprinted=0
debug=True
nPrint=12
mass_values=[] # create a vector where to store the computed values of mass
for e in LHEfile: # loop on the events
fourvectors_container=[] # to store the four-vector of this event
for p in e.particles: # loop on the particles of each event
if p.status == 1 and abs(p.id) == 11: # check it is a final state and is a e- or e+
lv=p.fourvector() # make four vector
if nprinted <nPrint:
if debug: print('[',nprinted,']',p.id, ':,',lv.px, lv.py, lv.pz, lv.e)
nprinted+=1
fourvectors_container.append(lv) # add the four-vector to the container
# end of loop on particles in the event
#nprinted=0
if len(fourvectors_container) == 4:
# start analizing the container of four vectors
_lv01 = fourvectors_container[0]+fourvectors_container[1]
_mass=_lv01.mass()
mass_values.append(_mass)
if nprinted <nPrint:
if debug: print(_lv01.px, _lv01.py, _lv01.pz, _lv01.e)
if debug: print(_mass)
nprinted+=1
# end loop on events
In [239]:
_mass
Out[239]:
In [245]:
hist, edges = np.histogram(mass_values,bins=arange(-1,120,1),normed=True)
p1 = figure(title="$m_{ee}$ Distribution Background",#tools="save",
background_fill_color="#E8DDCB")
p1.quad(top=hist, bottom=0, left=edges[:-1], right=edges[1:],
fill_color="#036564", line_color="#033649")
p1.xaxis.axis_label = 'm_(ee)'
p1.yaxis.axis_label = 'dσ/dm_{ee}'
show(p1)
In [248]:
plt.hist(mass_values,bins=arange(-1,120,1),density=True)
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{ee}$")
plt.ylabel("$d\sigma/dm_{ee}$")
plt.show()
Exercise
- repeat this plot for all the other possible pairs you can make with 4 electrons
- attempt to make a function that computes the mass of a specified pair $mass(i,j)$ which returns the mass of electron "i" and electron "j", for example mass(2,3) returns the mass of the 2nd and 3rd electron.
A third example: collect 4-plets as you go.¶
In [254]:
nprinted=0
debug=True
nPrint=12
mass_values=[] # create a vector where to store the computed values of mass
for e in LHEfile: # loop on the events
fourvectors_container=[] # to store the four-vector of this event
for p in e.particles: # loop on the particles of each event
if p.status == 1 and abs(p.id) == 11: # check it is a final state and is a e-
lv=p.fourvector() # make four vector
if nprinted <nPrint:
if debug: print('[',nprinted,']',p.id, ':,',lv.px, lv.py, lv.pz, lv.e)
nprinted+=1
fourvectors_container.append(lv) # add the four-vector to the container
# end of loop on particles in the event
#nprinted=0
if len(fourvectors_container) == 4:
# start analizing the container of four vectors
_lv01 = fourvectors_container[0]+fourvectors_container[1]+fourvectors_container[2]+fourvectors_container[3]
_mass=_lv01.mass()
mass_values.append(_mass)
if nprinted <nPrint:
if debug: print(_lv01.px, _lv01.py, _lv01.pz, _lv01.e)
if debug: print(_mass)
nprinted+=1
# end loop on events
In [252]:
hist, edges = np.histogram(mass_values,bins=arange(-1,1000,2),normed=True)
p1 = figure(title="$m_{eeee}$ Distribution Background",#tools="save",
background_fill_color="#E8DDCB")
p1.quad(top=hist, bottom=0, left=edges[:-1], right=edges[1:],
fill_color="#036564", line_color="#033649")
p1.xaxis.axis_label = 'm_(ee)'
p1.yaxis.axis_label = 'dσ/dm_{ee}'
show(p1)
In [255]:
plt.hist(mass_values,bins=arange(-1,1000,1),density=True)
plt.title("$m_{eeee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
In [257]:
plt.hist(mass_values,bins=arange(150,1000,1),density=True)
plt.title("$m_{eeee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
A 4th example: collect 4-plets as you go (and keep them for later).¶
In [276]:
nprinted=0
debug=True
nPrint=12
mass_values=[] # create a vector where to store the computed values of mass
for e in LHEfile: # loop on the events
fourvectors_container=[] # to store the four-vector of this event
for p in e.particles: # loop on the particles of each event
if p.status == 1 and abs(p.id) == 11: # check it is a final state and is a e-
lv=p.fourvector() # make four vector
if nprinted <nPrint:
if debug: print('[',nprinted,']',p.id, ':,',lv.px, lv.py, lv.pz, lv.e)
nprinted+=1
fourvectors_container.append(lv) # add the four-vector to the container
# end of loop on particles in the event
#nprinted=0
if len(fourvectors_container) == 4:
# start analizing the container of four vectors
_lv01 = fourvectors_container[0]+fourvectors_container[1]+fourvectors_container[2]+fourvectors_container[3]
_mass=_lv01.mass()
mass_values.append(_mass)
if nprinted <nPrint:
if debug: print(_lv01.px, _lv01.py, _lv01.pz, _lv01.e)
if debug: print(_mass)
nprinted+=1
# end loop on events
Keep the result and store it in a separate array for later plot
In [264]:
backgroud_mass=mass_values
run again for the signal and now keep the result and store it in a separate array for later plot.
In [271]:
signal_mass=mass_values
In [272]:
plt.hist([signal_mass,backgroud_mass],bins=arange(-1,1000,20),density=False)
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
In [186]:
plt.hist([signal_mass,backgroud_mass],bins=arange(-1,1000,20),density=False,weights=[1.15E-4*np.ones_like(signal_mass),1e-2*np.ones_like(backgroud_mass)])
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
In [273]:
plt.hist([signal_mass,backgroud_mass],bins=arange(-1,1000,20),log=True,density=False,weights=[1.15E-4*np.ones_like(signal_mass),1e-2*np.ones_like(backgroud_mass)])
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
In [274]:
plt.hist([signal_mass,backgroud_mass],bins=arange(-1,1000,2),log=True,stacked=True,density=False,weights=[1.15E-4*np.ones_like(signal_mass),1e-2*np.ones_like(backgroud_mass)])
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
In [275]:
LHEfile=lhef.readLHE("pp24l_nohiggs.lhe")
In [277]:
fullbackground_mass=mass_values
In [278]:
plt.hist([signal_mass,fullbackground_mass],bins=arange(-1,1000,2),log=True,stacked=True,density=False,weights=[1.15E-4*np.ones_like(signal_mass),6e-3*np.ones_like(backgroud_mass)])
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
In [279]:
plt.hist([signal_mass,fullbackground_mass],bins=arange(-1,1000,2),log=True,histtype='step',density=False,weights=[1.15E-4*np.ones_like(signal_mass),6e-3*np.ones_like(fullbackground_mass)])
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
In [201]:
plt.hist([signal_mass,fullbackground_mass],bins=arange(-1,1000,2),log=True,histtype='barstacked',density=False,weights=[1.15E-4*np.ones_like(signal_mass),16e-3*np.ones_like(fullbackground_mass)])
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
In [282]:
plt.hist([signal_mass,fullbackground_mass],bins=arange(-1,300,2),log=True,histtype='barstacked',density=False,weights=[1.15E-4*np.ones_like(signal_mass),16e-3*np.ones_like(fullbackground_mass)])
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
In [307]:
n, bins, patches = plt.hist([fullbackground_mass],bins=arange(-1,1000,2),log=True,histtype='step',density=False,weights=[1.0E-2*np.ones_like(fullbackground_mass)])
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
In [304]:
print(bins[60:66],'--->',n[60:66])
In [300]:
n, bins, patches = plt.hist([signal_mass],bins=arange(-1,1000,2),log=True,histtype='step',density=False,weights=[1.15E-4*np.ones_like(signal_mass)])
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
In [302]:
print(bins[60:66],'--->',n[60:66])
In [281]:
plt.hist([fullbackground_mass],bins=arange(-1,110,2),log=True,histtype='step',density=False,weights=[1.15E-4*np.ones_like(fullbackground_mass)])
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{eeee}$")
plt.ylabel("$d\sigma/dm_{eeee}$")
plt.show()
A second example: collect pairs as you go.¶
In [203]:
nprinted=0
debug=True
nPrint=12
mass_values=[] # create a vector where to store the computed values of mass
for e in LHEfile: # loop on the events
fourvectors_container=[] # to store the four-vector of this event
for p in e.particles: # loop on the particles of each event
if p.status == 1 and abs(p.id) == 11: # check it is a final state and is a e-
lv=p.fourvector() # make four vector
if nprinted <nPrint:
if debug: print('[',nprinted,']',p.id, ':,',lv.px, lv.py, lv.pz, lv.e)
nprinted+=1
fourvectors_container.append(lv) # add the four-vector to the container
# end of loop on particles in the event
#nprinted=0
if len(fourvectors_container) == 4:
# start analizing the container of four vectors
_lv01 = fourvectors_container[0]+fourvectors_container[1]
_mass=_lv01.mass()
mass_values.append(_mass)
if nprinted <nPrint:
if debug: print(_lv01.px, _lv01.py, _lv01.pz, _lv01.e)
if debug: print(_mass)
nprinted+=1
# end loop on events
In [204]:
hist, edges = np.histogram(mass_values,bins=arange(-1,120,1),normed=True)
p1 = figure(title="$m_{ee}$ Distribution Background",#tools="save",
background_fill_color="#E8DDCB")
p1.quad(top=hist, bottom=0, left=edges[:-1], right=edges[1:],
fill_color="#036564", line_color="#033649")
p1.xaxis.axis_label = 'm_(ee)'
p1.yaxis.axis_label = 'dσ/dm_{ee}'
show(p1)
In [136]:
plt.hist(mass_values,bins=arange(-1,120,1),density=True)
plt.title("$m_{ee}$ Distribution Signal")
plt.xlabel("$m_{ee}$")
plt.ylabel("$d\sigma/dm_{ee}$")
plt.show()