forked from dally96/DisplacedTauAnalysis
-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathdisEleReco.py
More file actions
124 lines (99 loc) · 5.52 KB
/
disEleReco.py
File metadata and controls
124 lines (99 loc) · 5.52 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
import uproot
import scipy
import matplotlib as mpl
import awkward as ak
import numpy as np
import math
import ROOT
import array
import pandas as pd
stauFile = uproot.open("Staus_M_100_100mm_13p6TeV_Run3Summer22_NanoAOD.root")
Events = stauFile["Events"]
eleBranches = {}
genPartBranches = {}
can = ROOT.TCanvas("can", "can")
dxy_bin = [-30., -25., -20., -15., -10., -5., -4.5, -4., -3.5, -3., -2.5, -2., -1.5, -1., -0.5, 0,
0.5, 1., 1.5, 2., 2.5, 3., 3.5, 4., 4.5, 5., 10., 15., 20., 25., 30.]
dxy_bin = array.array('d', dxy_bin)
for key in Events.keys():
if "Electron_" in key:
eleBranches[key] = Events[key].array()
if "GenPart_" in key:
genPartBranches[key] = Events[key].array()
def isMotherStau(evt, part_idx):
#Find its mother
mother_idx = genPartBranches["GenPart_genPartIdxMother"][evt][part_idx]
#If it's a stau, we're done
if (abs(genPartBranches["GenPart_pdgId"][evt][mother_idx]) == 1000015):
return 1
#If it's a tau, go further up the chain
if (abs(genPartBranches["GenPart_pdgId"][evt][mother_idx]) == 15):
return isMotherStau(evt, mother_idx)
#Else, it's not what we're looking for
else:
return 0
def isMotherStauIdx(evt, part_idx):
mother_idx = genPartBranches["GenPart_genPartIdxMother"][evt][part_idx]
if (abs(genPartBranches["GenPart_pdgId"][evt][mother_idx]) == 1000015):
return mother_idx
else:
return isMotherStauIdx(evt, mother_idx)
# This chooses the gen particles which are electrons which are decayed from displaced taus
genDisEle = []
genStauIdx = []
for evt in range(len(genPartBranches["GenPart_pdgId"])):
genDisEle_evt = []
genStauIdx_evt = []
if (len(genPartBranches["GenPart_pdgId"][evt]) != 0):
for part_idx in range(len(genPartBranches["GenPart_pdgId"][evt])):
if (abs(genPartBranches["GenPart_pdgId"][evt][part_idx]) == 11):
if (isMotherStau(evt, part_idx)):
genDisEle_evt.append(part_idx)
genStauIdx_evt.append(isMotherStauIdx(evt, part_idx))
genDisEle.append(genDisEle_evt)
genStauIdx.append(genStauIdx_evt)
genDisElePt = genPartBranches["GenPart_pt"][genDisEle]
genDisEleEta = genPartBranches["GenPart_eta"][genDisEle]
genDisElePhi = genPartBranches["GenPart_phi"][genDisEle]
def dxy(evt, muon_idx, stau_idx):
return (genPartBranches["GenPart_vertexY"][evt][muon_idx] - genPartBranches["GenPart_vertexY"][evt][stau_idx]) * np.cos(genPartBranches["GenPart_phi"][evt][muon_idx]) - (genPartBranches["GenPart_vertexX"][evt][muon_idx] - genPartBranches["GenPart_vertexX"][evt][stau_idx]) * np.sin(genPartBranches["GenPart_phi"][evt][muon_idx])
genPartBranches["GenEle_dxy"] = []
for evt in range(len(genPartBranches["GenPart_pdgId"])):
genDisEledxy_evt = np.zeros(len(genPartBranches["GenPart_pdgId"][evt]))
if (len(genPartBranches["GenPart_pdgId"][evt]) != 0):
for idx in range(len(genDisEle[evt])):
dxy_val = dxy(evt, genDisEle[evt][idx], genStauIdx[evt][idx])
genDisEledxy_evt[genDisEle[evt][idx]] = dxy_val
genPartBranches["GenEle_dxy"].append(ak.Array(genDisEledxy_evt))
genPartBranches["GenEle_dxy"] = ak.Array(genPartBranches["GenEle_dxy"])
#print(genPartBranches["GenMuon_dxy"])
genDisEledxy = genPartBranches["GenEle_dxy"][genDisEle]
genDisElePt = genDisElePt[genDisElePt > 20]
genDisEleEta = genDisEleEta[genDisElePt > 20]
genDisEledxy = genDisEledxy[genDisElePt > 20]
filteredRecoEle = [np.array([val for val in subarr if val in values]) for subarr, values in zip(eleBranches["Electron_genPartIdx"], genDisEle)]
filteredRecoEleIndices = [np.where(np.isin(subarr, values))[0] for subarr, values in zip(eleBranches["Electron_genPartIdx"], genDisEle)]
filteredRecoElePt = genPartBranches["GenPart_pt"][filteredRecoEle]
filteredRecoEleEta = genPartBranches["GenPart_eta"][filteredRecoEle]
filteredRecoEledxy = genPartBranches["GenEle_dxy"][filteredRecoEle]
filteredRecoElePt = filteredRecoElePt[filteredRecoElePt > 20]
filteredRecoEleEta = filteredRecoEleEta[filteredRecoElePt > 20]
filteredRecoEledxy = filteredRecoEledxy[filteredRecoElePt > 20]
disEle_pteff = ROOT.TEfficiency("disEle_pteff", "Staus_M_100_100mm_13p6TeV_Run3Summer22; pt [GeV]; Fraction of gen dis ele which are recon'd", 25, 0, 100)
for i in range(25):
disEle_pteff.SetTotalEvents(i + 1, len(ak.flatten(genDisElePt[(genDisElePt > (i * 4)) & (genDisElePt < ((i + 1) * 4))])))
disEle_pteff.SetPassedEvents(i + 1, len(ak.flatten(filteredRecoElePt[(filteredRecoElePt > (i * 4)) & (filteredRecoElePt < ((i + 1) * 4))])))
disEle_etaeff = ROOT.TEfficiency("disEle_etaeff", "Staus_M_100_100mm_13p6TeV_Run3Summer22; #eta; Fraction of gen dis ele which are recon'd", 16, -2.4, 2.4)
for i in range(16):
disEle_etaeff.SetTotalEvents(i + 1, len(ak.flatten(genDisEleEta[(genDisEleEta > (-2.4 + 0.3 * i)) & (genDisEleEta < (-2.4 + (i + 1) * 0.3))])))
disEle_etaeff.SetPassedEvents(i + 1, len(ak.flatten(filteredRecoEleEta[(filteredRecoEleEta > (i * 0.3 - 2.4)) & (filteredRecoEleEta < ((i + 1) * 0.3 - 2.4))])))
disEle_dxyeff = ROOT.TEfficiency("disEle_dxyeff", "Staus_M_100_100mm_13p6TeV_Run3Summer22; d_{xy} [cm]; Fraction of gen dis ele which are recon'd", len(dxy_bin) - 1, dxy_bin)
for i in range(len(dxy_bin) - 1):
disEle_dxyeff.SetTotalEvents(i + 1, len(ak.flatten(genDisEledxy[(genDisEledxy >= dxy_bin[i]) & (genDisEledxy < dxy_bin[i + 1])])))
disEle_dxyeff.SetPassedEvents(i + 1, len(ak.flatten(filteredRecoEledxy[(filteredRecoEledxy >= dxy_bin[i]) & (filteredRecoEledxy < dxy_bin[i + 1])])))
disEle_pteff.Draw()
can.SaveAs("disEle_recoPtEff.pdf")
disEle_etaeff.Draw()
can.SaveAs("disEle_recoEtaEff.pdf")
disEle_dxyeff.Draw()
can.SaveAs("disEle_recoDxyEff.pdf")