CERN Z Boson Identification & Understanding
INTRODUCTION
In the following article, I shall be analysing data received from the CMS detectors at CERN. Its purpose is to detect the different particles that may appear after a collision of a group of protons in the Large Hadron Collider; this allows the scientist to study the constitution of our universe, dark matter, et cetera. The data collected used in our module considers purely muons originating from the collided particles. Muons are elementary particles classified as leptons. Thus it is believed to have no sub-structure. With the data we collected, we shall be primarily extracting information from each muon to get the invariant mass of the larger particles. The invariant mass is the mass that a body has once at rest; the particles measured are going at high velocities. Therefore, we do not have their ‘real mass’ as it is distorted aka the relativistic mass. The goal of the data analysis is to analyse individual, more significant particles, in our case, we shall be looking at the z boson or the z zero particle. CERN discovered this boson in the early 1980s. Much like the W boson, the z zero particle carries the weak force, unlike his counterpart(W±) this boson is neutral.
HISTORY
To further understand the particle that we shall be converging on, we shall take a look at the history of it. At CERN, scientists realised that by emitting a W boson, the weak force could change; for example, the proton’s charge by altering its quark configuration. Sidney Bludman in 1958 suggested that there could be another actor in the weak force, but this time a neutral counterpart. Once experiments initiated, scientist shot neutrinos into a bubble chamber in which the neutrinos started interacting, thus proving the there was a weak neutral current: the Z boson. Furthermore, after the creation of the super proton synchrotron (end of the 1970s), the proton-antiproton collider they were able to prove the presence of both the W and Z bosons empirically, thanks to the UA1 which was a particle detector in the CERN SPS.
The Data
The data received from the CMS comes in a root file which contains information about over 57 million collisions. These events produced muons which have been documented; the datafile lets us know their charge, transverse momentum, phi and theta angle once collided. A muon is either positively or negatively charged (1 or -1). The transverse momentum of the muon equates to the projection of the momentum(3d vector) on the transverse plane, making a 2d figure. The theta angle of the muon is equal to the angle between the y-axis and the particle projection, whereas the phi angle is equal to the angle between the particle projection and the x-axis. All angles were measured in radians.
%matplotlib inline
import numpy as np
import matplotlib.mlab as mlab
import matplotlib.pyplot as plt
import ROOT
import math
import datetime
read = ROOT.TFile.Open("Run2012BC_DoubleMuParked.root")
parse = read.Get("Events")
phi={}
theta={}
print "start"
df=[]
start = datetime.datetime.now()
k=0for i in parse:
k+=1
if k>10000:break
if k%1000000 == 0:
print(k)
for g in range(i.nMuon):
fig1 = int(-(i.Muon_phi[g]-np.pi)/(np.pi/100))
cfig1= int(((fig1 * np.pi/100)- np.pi + (np.pi/200))*100)
fig2 = int(-(i.Muon_eta[g]-np.pi)/(np.pi/100))
cfig2 = int(((fig2 * np.pi/100)- np.pi + (np.pi/200))*100)
if cfig1 in phi:
phi[cfig1]+=1
else:
phi[cfig1]=1
if fig2 in theta:
theta[fig2]+=1
else:
theta[fig2]=1
print(datetime.datetime.now()-start)
f=plt.gcf()
plt.title("Phi")
f.set_size_inches(9.25, 5.25)
plt.bar([i for i in phi],[i for i in phi.values()], color="blue")
plt.show()
g=plt.figure(2)
plt.title("Theta")
plt.bar([i for i in theta], [i for i in theta.values()], color="red")
g.show()
plt.xlabel('Values')
plt.ylabel('Frequency')
plt.show()
Invariant Mass
The invariant mass plot gives us a very strong indication of the particle that we are dealing with. We notice that at a particular time, the invariant mass rises majorly, which could indicate the presence of a particle. In our dataset, we notice that there are various peaks, but we shall only focus on one in this notebook, namely the Z boson (as I mentioned previously) which is the last peak of the invariant mass graph.
def massm(pt,et,ph,ms):
mass1 = ms[0]
px1 = vector(pt[0],et[0],ph[0])[0]
py1 = vector(pt[0],et[0],ph[0])[1]
pz1 = vector(pt[0],et[0],ph[0])[2]
energy1 = math.sqrt(mass1**2 + px1**2 + py1**2 + pz1**2)
mass2 = ms[1]
px2 = vector(pt[1],et[1],ph[1])[0]
py2 = vector(pt[1],et[1],ph[1])[1]
pz2 = vector(pt[1],et[1],ph[1])[2]
energy2 = math.sqrt(mass2**2 + px2**2 + py2**2 + pz2**2)
mass1 = energy1**2 + 2*energy1*energy2 + energy2**2
mass2 = (px1+px2)**2 + (py1+py2)**2 + (pz1+pz2)**2
mass = math.sqrt(mass1-mass2)
return massphi=[]
theta=[]
pt=[]
mass=[]
%pylab inline
mass2=[]
print "start"
k=0
for i in parse:
k+=1
if k>100000:break
if k%1000000 == 0:
print(k)
for g in range(i.nMuon):
phi.append(i.Muon_phi[g])
theta.append(i.Muon_eta[g])
tt = i.Muon_pt[g]
if tt<150: pt.append(tt)
if i.nMuon!=2: continue
if i.Muon_charge[0] * i.Muon_charge[1] < 0:
qv0 = ROOT.Math.PtEtaPhiMVector(i.Muon_pt[0], i.Muon_eta[0], i.Muon_phi[0], i.Muon_mass[0])
qv1 = ROOT.Math.PtEtaPhiMVector(i.Muon_pt[1], i.Muon_eta[1], i.Muon_phi[1], i.Muon_mass[1])
ms = (qv0 + qv1).mass()
if ms<20:
mass2.append(ms)
if ms<120:
mass.append(ms)f=plt.figure(1)
plt.title("Phi")
plt.hist(phi, color="blue")
f.set_size_inches(9.25, 5.25)
f.show()
g=plt.figure(2)
plt.title("Theta")
plt.hist(theta, color="red")
g.set_size_inches(9.25, 5.25)
g.show()
e=plt.figure(3)
plt.title("Pt")
plt.hist(pt, color="green", bins=10)
e.set_size_inches(9.25, 5.25)
e.show()
t=plt.figure(4)
plt.title("invariant mass")
plt.hist(mass, color="yellow", bins=150)
plt.yscale("log")
t.set_size_inches(9.25, 5.25)
t.show()
h=plt.figure(5)
plt.title("invariant mass")
plt.hist(mass2, color="yellow", bins='auto')
plt.yscale("log")
h.set_size_inches(9.25, 5.25)
h.show()
mass=[]
bmass = []
gmass=[]
k=0
g=0
for i in parse:
k+=1
if k%100000 == 0:
print(k)
if k>1000000:break
if i.nMuon!=2: continue
qv0 = ROOT.Math.PtEtaPhiMVector(i.Muon_pt[0], i.Muon_eta[0], i.Muon_phi[0], i.Muon_mass[0])
qv1 = ROOT.Math.PtEtaPhiMVector(i.Muon_pt[1], i.Muon_eta[1], i.Muon_phi[1], i.Muon_mass[1])
ms = (qv0 + qv1).mass()
p = math.sqrt((qv0.X()+qv1.X())**2 + (qv0.Y()+qv1.Y())**2 + (qv0.Z()+qv1.Z())**2)
if p > 10:
if i.Muon_charge[0] * i.Muon_charge[1] > 0:
if ms<120:
bmass.append(ms)
elif i.Muon_charge[0] * i.Muon_charge[1] < 0:
if ms<120:
gmass.append(ms)
def mean_clean(cleanhist,xhist,start,end):
xv = xhist[1].tolist()
yv = cleanhist
k=0
sum = 0
for i,j in zip(xv,yv):
if i < start: continue
if i > end: continue
i = (xv[xv.index(i) + 1] + i)/2
sum += i * j
k += j
return sum/kBackground
The raw data/graph, however, can be relatively misleading. In the CMS, the conditions are close to perfect, but not exactly there yet, meaning that there are different possibilities of disturbance that can affect the detection, or the quality of the data receives. In our case for example, even if we have a vast amount of data the disturbances are still noticeable shown by the background in the first histogram of the invariant mass because data is a collection of many collisions and Z zero is not always present. To remove a reasonable amount of background, we removed 2 significant sources of background. Our first variable was the mass of the particle: if the invariant mass of the possible particle was not noteworthy, we decided that it was probably not a particle but possibly an offshoot of a prior collision. The second variable that we considered was the momentum. We concluded that if the momentum of a possible particle was too low, then it wasn’t necessary for it to be included in the calculations (since the Z zero has a relatively large mass, I would expect the momentum to be significant).
b=[]
for i in bmass:
b.append(2.4)
h = plt.figure(1)
plt.yscale("log")
one = plt.hist(gmass, bins=500, facecolor='green')
two = plt.hist(bmass,bins=500, facecolor='yellow', weights=b, alpha=0.7)
plt.axis([60, 120, 0, 3000])
h.set_size_inches(9.25, 5.25)
plt.show()
The Subtraction
Once we are done estimating/identifying the background, we have to find a way to remove from the raw histogram. To achieve this, we try to align the background’s histogram as much as possible to the raw histogram to see how much the background is affecting our data. After the alignment is done, we construct a bar chart by subtracting the arrays of the y coordinates of both histograms (raw — background). Once this is done, the peak where the Z boson is supposed to be should ideally be alone and parallel to the x-axis.
Mass of Z zero
We want to calculate the Z boson’s mass. To achieve this, we will first need to remove the background(explained in the paragraph above). This eradication of erroneous sources allows the peak to be parallel with the x-axis thus easier to calculate. The removal of the disturbances shows us a clean figure of the boson detection. Its mass is equal to the x coordinate of the vertex of this peak. To retrieve this information, we calculate the mean value of the bins included inside the peak, which should give us an approximately correct value of its mass in GeV/c². If we haven’t committed any considerable mistakes, the value that we get should be close to the legitimate value of the Z boson, which is around 91.18 GeV/c². Of course, if we consider more possible erroneous sources, other than the low momentum and the insignificant mass, the value which we calculate will continuously move closer to the real one.
FWHM
After having calculated the mass of the particle, another interesting measurement would be the full width at half maximum(FWHM) value. This value represents the decay width of a particle which will calculate by considering the remaining background as 0 and calculating the distance between the halfway points. The distance represents the decay width which is directly related to the decay lifetime as well as being the direct inverse of the particle’s resonance. The value that we receive seems to be rather large compared to its supposed value, which is 2.49 GeV/c². The environment of the detector causes this divergence. The particle, before being detected, goes through multiple layers of substance/metals, so its ‘flight path’ diverges. Thus the particle’s properties (e.g. velocity, direction, et cetera) can be altered. These minimal litigating circumstances or multiple scatterings will, therefore, broaden the results perceived by the detector. Ipso facto the peaks of the invariant masses’ histograms become more massive and may skew the results.
bar = []
for d in range(len(one[1])):
if d==0: continue
middle = (one[1][d]+one[1][d-1])/2
bar.append(middle)
new=[]
for i,j in zip(one[0],two[0]):
new.append(i-(j*2.3))
point = mean_clean(new,one,87,94)
print(point)
big=0
k=0
for i in new:
k+=1
if k<80:continue
if i>big:
big=i
middle = big/2
num=0
k=0
high = new.index(big)
for i in new:
k+=1
if k<high:continue
if abs(middle-i)<abs(middle-num):
num=i
first = one[1][new.index(num)]
num=0
k=0
for i in new:
k+=1
if k<70:continue
if k>high:continue
if abs(middle-i)<abs(middle-num):
num=i
second = one[1][new.index(num)]
dif = first - second
print dif
p3=[first,second]
p4=[middle, middle]
p1=[0,100000]
p2=[point,point]
diff=plt.bar(bar, height=(one[0]-two[0]))
plt.plot(p2,p1,color="red")
plt.plot(p3,p4,color="red")
plt.axis([70, 110, 0, 4000])
plt.show()
