In [1]:
import ROOT
%jsroot on
import os
import glob
import uproot
import numpy as np
import matplotlib.pyplot as plt
In [2]:
def format_num(n):
if n > 1e7: return f'{n/1e6:.0f}M'
if n > 1e6: return f'{n/1e6:.1f}M'
if n > 1e4: return f'{n/1e3:.0f}k'
if n > 1e3: return f'{n/1e3:.1f}k'
else: return f'{n}'
# format_num(12345678)
Plans for next week¶
- compare phi distribution with MC
total and 1 uniform and 1 non-uniform runs - how well are they reproduced - QA for 244456 (and 244453)
- check triggers settings - it should not be the case as the pT spectra does not differ
- check global event properties: event multiplicities, event vertex distr.
- check RCT & logbook
- rather small stats of pp@5.02TeV -- what with the reference for PbPb?
check approaches in PbPb@5.02TeV papers
Compare $\phi$ distribution with MC¶
LHC15n only¶
In [10]:
f = ROOT.TFile('../../HF-jets/ana_results/iter2/LHC15n/AnalysisResults.root')
t = ROOT.gROOT.FindObject('JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_allJets')
In [11]:
# ROOT.gROOT.Reset()
c1 = ROOT.TCanvas( 'c1', 'Example', 200, 10, 700, 500 )
c1.Draw()
h_jet_phi = ROOT.TH1D('h_jet_phi', 'h_jet_phi', 100, 0, 6.283)
h_jet_phi.GetXaxis().SetTitle('jet phi')
t.Draw('Jet_Phi >> h_jet_phi', '', 'e1x0')
# h_jet_phi.DrawCopy()
Out[11]:
In [12]:
ROOT.gROOT.Reset()
c1_2 = ROOT.TCanvas( 'c1_2', '', 200, 10, 700, 500 )
c1_2.Draw()
h_track_phi = ROOT.TH1D('h_track_phi', 'h_track_phi', 100, 0, 6.283)
h_track_phi.GetXaxis().SetTitle('track phi')
t.Draw('Jet_Track_Phi >> h_track_phi', '', 'e1x0')
# h_phi.GetYaxis().SetRangeUser(6000,8000)
# c1.Draw()
Out[12]:
Is phi oscillations amplitude dependend on $p_T$ ???¶
In [17]:
nbins = 90
f = ROOT.TFile(f'../../HF-jets/ana_results/iter2/LHC16h3/ptbin15/AnalysisResults.root')
t = ROOT.gROOT.FindObject('JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_udsgJets')
c1 = ROOT.TCanvas( 'c1', 'Example', 900, 500 )
c1.Draw()
h1 = ROOT.TH1D('h1', 'h1', nbins, 0, 6.283)
h2 = ROOT.TH1D('h2', 'h2', nbins, 0, 6.283)
h3 = ROOT.TH1D('h3', 'h3', nbins, 0, 6.283)
h1.SetLineColor(ROOT.kBlue)
h2.SetLineColor(ROOT.kRed)
h3.SetLineColor(ROOT.kBlack)
t.Draw(f'Jet_Track_Phi >> h1', 'Jet_Track_Pt < 1', 'goff')
t.Draw(f'Jet_Track_Phi >> h2', 'Jet_Track_Pt > 1 & Jet_Track_Pt < 3', 'goff')
t.Draw(f'Jet_Track_Phi >> h3', 'Jet_Track_Pt > 3', 'goff')
print(f'N entries: {format_num(h1.GetEntries())}, {format_num(h2.GetEntries())}, {format_num(h3.GetEntries())}')
h1.SetTitle('pT < 1')
h2.SetTitle('1 < pT < 3')
h3.SetTitle('pT > 3')
h1.Scale(1/h1.Integral())
h2.Scale(1/h2.Integral())
h3.Scale(1/h3.Integral())
h3.Draw('le')
h2.Draw('le,same')
h1.Draw('le,same')
c1.GetPad(0).BuildLegend(0.55, 0.775, 0.78, 0.935)
Out[17]:
whole LHC15n vs MC (LHC16h3)¶
No cut¶
In [18]:
nbins = 90
cut = ''
f_data = ROOT.TFile('../../HF-jets/ana_results/iter2/LHC15n/AnalysisResults.root')
t_data = ROOT.gROOT.FindObject('JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_allJets')
c1 = ROOT.TCanvas( 'c1', 'Example', 900, 500 )
c1.Draw()
h_track_phi = ROOT.TH1D('h_data', 'h_data', nbins, 0, 6.283)
h_track_phi.SetTitle('data LHC15n')
h_track_phi.GetXaxis().SetTitle('track phi')
h_track_phi.SetLineColor(2)
h_track_phi.SetLineWidth(4)
t_data.Draw('Jet_Track_Phi >> h_data', cut, 'e,goff')
h_track_phi.Scale(1./h_track_phi.Integral())
h_track_phi.SetMinimum(1/nbins*0.8)
h_track_phi.SetMaximum(1/nbins*1.2)
h_track_phi.Draw('le')
ROOT.gStyle.SetPalette(55)
for i,ptbin in enumerate([ '10', '15']):
f_mc = ROOT.TFile(f'../../HF-jets/ana_results/iter2/LHC16h3/ptbin{ptbin}/AnalysisResults.root')
t_mc = ROOT.gROOT.FindObject('JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_udsgJets')
h_track_phi_mc = ROOT.TH1D(f'h_mc_ptbin{ptbin}', f'h_mc_ptbin{ptbin}', nbins, 0, 6.283)
h_track_phi_mc.SetLineColor(ROOT.gStyle.GetColorPalette(i*150))
h_track_phi_mc.SetLineWidth(2)
t_mc.Draw(f'Jet_Track_Phi >> h_mc_ptbin{ptbin}', cut, 'e,goff')
h_track_phi_mc.Scale(1./h_track_phi_mc.Integral())
h_track_phi_mc.DrawCopy('le,same')
c1.GetPad(0).BuildLegend(0.55, 0.775, 0.78, 0.935)
Out[18]:
Low pT¶
In [19]:
nbins = 90
cut = 'Jet_Track_Pt < 1'
f_data = ROOT.TFile('../../HF-jets/ana_results/iter2/LHC15n/AnalysisResults.root')
t_data = ROOT.gROOT.FindObject('JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_allJets')
c1 = ROOT.TCanvas( 'c1', 'Example', 900, 500 )
c1.Draw()
h_track_phi = ROOT.TH1D('h_data', 'h_data', nbins, 0, 6.283)
h_track_phi.SetTitle('data LHC15n')
h_track_phi.GetXaxis().SetTitle('track phi')
h_track_phi.SetLineColor(2)
h_track_phi.SetLineWidth(4)
t_data.Draw('Jet_Track_Phi >> h_data', cut, 'e,goff')
h_track_phi.Scale(1./h_track_phi.Integral())
h_track_phi.SetMinimum(1/nbins*0.8)
h_track_phi.SetMaximum(1/nbins*1.2)
h_track_phi.Draw('le')
ROOT.gStyle.SetPalette(55)
for i,ptbin in enumerate([ '10', '15']):
f_mc = ROOT.TFile(f'../../HF-jets/ana_results/iter2/LHC16h3/ptbin{ptbin}/AnalysisResults.root')
t_mc = ROOT.gROOT.FindObject('JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_udsgJets')
h_track_phi_mc = ROOT.TH1D(f'h_mc_ptbin{ptbin}', f'h_mc_ptbin{ptbin}', nbins, 0, 6.283)
h_track_phi_mc.SetLineColor(ROOT.gStyle.GetColorPalette(i*150))
h_track_phi_mc.SetLineWidth(2)
t_mc.Draw(f'Jet_Track_Phi >> h_mc_ptbin{ptbin}', cut, 'e,goff')
h_track_phi_mc.Scale(1./h_track_phi_mc.Integral())
h_track_phi_mc.DrawCopy('le,same')
c1.GetPad(0).BuildLegend(0.55, 0.775, 0.78, 0.935)
Out[19]:
Per run¶
No cut¶
In [20]:
nbins = 90
cut = ''
for run in [r[3:] for r in os.listdir('../../HF-jets/ana_results/iter2/LHC15n/myOutputDir/') if r.startswith('000')]:
print(f'run = {run}')
f_data = ROOT.TFile(f'../../HF-jets/ana_results/iter2/LHC15n/myOutputDir/000{run}/AnalysisResults.root')
t_data = ROOT.gROOT.FindObject('JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_allJets')
# ROOT.gROOT.Reset()
c1 = ROOT.TCanvas( 'c1'+run, 'c1'+run, 900, 500 )
c1.Draw()
h_track_phi = ROOT.TH1D('h_data', 'h_data', nbins, 0, 6.283)
h_track_phi.GetXaxis().SetTitle(f'track phi')
h_track_phi.SetTitle(f'track phi {run}')
h_track_phi.SetLineColor(2)
h_track_phi.SetLineWidth(3)
t_data.Draw('Jet_Track_Phi >> h_data', cut, 'e,goff')
n = h_track_phi.GetEntries()
print(f'\t N entries = {format_num(n)} \t err ~ {np.sqrt(n/nbins)/(n/nbins)*100:.1f}%')
h_track_phi.Scale(nbins/h_track_phi.Integral())
h_track_phi.SetMinimum(nbins/nbins*0.91)
h_track_phi.SetMaximum(nbins/nbins*1.09)
h_track_phi.DrawCopy('le')
ROOT.gStyle.SetPalette(55)
for i,ptbin in enumerate([ '10', '15']):
try:
f_mc = ROOT.TFile(f'../../HF-jets/ana_results/iter2/LHC16h3/ptbin{ptbin}/myOutputDir/{run}/AnalysisResults.root')
t_mc = ROOT.gROOT.FindObject('JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_udsgJets')
h_track_phi_mc = ROOT.TH1D(f'h_mc_ptbin{ptbin}', f'h_mc_ptbin{ptbin}', nbins, 0, 6.283)
h_track_phi_mc.SetLineColor(ROOT.gStyle.GetColorPalette(i*200))
h_track_phi.SetLineWidth(2)
t_mc.Draw(f'Jet_Track_Phi >> h_mc_ptbin{ptbin}', cut, 'e,goff')
except:
f_mc.Close()
continue
n = float(h_track_phi_mc.GetEntries())
print(f'\t N entries ({ptbin}) = {format_num(n)} \t err ~ {np.sqrt(n/nbins)/(n/nbins)*100:.1f}%')
h_track_phi_mc.Scale(nbins/h_track_phi_mc.Integral())
opt='le,same'
h_track_phi_mc.DrawCopy(opt)
f_mc.Close()
c1.GetPad(0).BuildLegend(0.55, 0.775, 0.78, 0.935)
f_data.Close()
Low pT¶
In [21]:
nbins = 90
cut = 'Jet_Track_Pt < 1'
for run in [r[3:] for r in os.listdir('../../HF-jets/ana_results/iter2/LHC15n/myOutputDir/') if r.startswith('000')]:
print(f'run = {run}')
f_data = ROOT.TFile(f'../../HF-jets/ana_results/iter2/LHC15n/myOutputDir/000{run}/AnalysisResults.root')
t_data = ROOT.gROOT.FindObject('JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_allJets')
# ROOT.gROOT.Reset()
c1 = ROOT.TCanvas( 'c1'+run, 'c1'+run, 900, 500 )
c1.Draw()
h_track_phi = ROOT.TH1D('h_data', 'h_data', nbins, 0, 6.283)
h_track_phi.GetXaxis().SetTitle(f'track phi')
h_track_phi.SetTitle(f'track phi {run}')
h_track_phi.SetLineColor(2)
h_track_phi.SetLineWidth(3)
t_data.Draw('Jet_Track_Phi >> h_data', cut, 'e,goff')
n = h_track_phi.GetEntries()
print(f'\t N entries = {format_num(n)} \t err ~ {np.sqrt(n/nbins)/(n/nbins)*100:.1f}%')
h_track_phi.Scale(nbins/h_track_phi.Integral())
h_track_phi.SetMinimum(nbins/nbins*0.91)
h_track_phi.SetMaximum(nbins/nbins*1.09)
h_track_phi.DrawCopy('le')
ROOT.gStyle.SetPalette(55)
for i,ptbin in enumerate([ '10', '15']):
try:
f_mc = ROOT.TFile(f'../../HF-jets/ana_results/iter2/LHC16h3/ptbin{ptbin}/myOutputDir/{run}/AnalysisResults.root')
t_mc = ROOT.gROOT.FindObject('JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_udsgJets')
h_track_phi_mc = ROOT.TH1D(f'h_mc_ptbin{ptbin}', f'h_mc_ptbin{ptbin}', nbins, 0, 6.283)
h_track_phi_mc.SetLineColor(ROOT.gStyle.GetColorPalette(i*200))
h_track_phi.SetLineWidth(2)
t_mc.Draw(f'Jet_Track_Phi >> h_mc_ptbin{ptbin}', cut, 'e,goff')
except:
f_mc.Close()
continue
n = float(h_track_phi_mc.GetEntries())
print(f'\t N entries ({ptbin}) = {format_num(n)} \t err ~ {np.sqrt(n/nbins)/(n/nbins)*100:.1f}%')
h_track_phi_mc.Scale(nbins/h_track_phi_mc.Integral())
opt='le,same'
h_track_phi_mc.DrawCopy(opt)
f_mc.Close()
c1.GetPad(0).BuildLegend(0.55, 0.775, 0.78, 0.935)
f_data.Close()
QA for 244456 and 244453 (Ad 2)¶
Event properties¶
In [3]:
%jsroot off
In [4]:
data_files = glob.glob('../../HF-jets/ana_results/iter2/LHC15n/myOutputDir/*/AnalysisResults.root')
for i,run_file in [f for f in enumerate(data_files)]:
# ROOT.gROOT.Reset()
print(f'{run_file}, {i+1}/{len(data_files)}')
run = run_file.split('/')[-2]
f_run = ROOT.TFile(run_file)
f_run.Close()
# kernel sometimes crashes...
# continue
f_run = ROOT.TFile(run_file)
# continue
t_run = ROOT.gROOT.FindObject('JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_allJets')
c = ROOT.TCanvas( 'c_'+run, 'c_'+run, 200, 10, 950, 400 )
c.Draw()
pad1 = ROOT.TPad("pad1"+run, "pad1"+run, 0, 0, 0.5, 1)
# pad1.SetGridx()
pad1.Draw()
# Lower ratio plot is pad2
c.cd() # returns to main canvas before defining pad2
pad2 = ROOT.TPad("pad2"+run, "pad2"+run, 0.5, 0, 1, 1)
# pad2.SetBottomMargin(0.2)
# pad2.SetGridx()
pad2.Draw()
cut = 'Jet_Pt > 5'
pad1.cd()
h_2d = ROOT.TH2D('h_2d', 'h_2d', 100, 0.29, 0.38, 100, 0.02, 0.13)
h_2d.SetTitle(run)
h_2d.GetXaxis().SetTitle('vertex X')
h_2d.GetYaxis().SetTitle('vertex Y')
t_run.Draw('Event_Vertex_X:Event_Vertex_Y >> h_2d', cut, 'goff')
h_2d.DrawCopy('colz')
# c2 = ROOT.TCanvas( 'c_2_'+run, 'c_2_'+run, 200, 10, 700, 500 )
# c2.Draw()
pad2.cd()
pad2.SetLogy(1)
# h_jet_pt = ROOT.TH1D('h_jet_pt', 'h_jet_pt', 100, 0, 100)
# h_jet_pt.SetTitle(run)
# h_jet_pt.GetXaxis().SetTitle('jet pt')
# t_run.Draw('Jet_Pt >> h_jet_pt', cut, 'goff')
# h_jet_pt.DrawCopy()
h1 = ROOT.TH1D('h1', 'h1', 100, 0, 100)
h1.SetTitle(run)
h1.GetXaxis().SetTitle('event multiplicity')
t_run.Draw('Event_Multiplicity >> h1', cut, 'goff')
h1.DrawCopy()
print('done')
f_run.Close()
# canvases.append(c)
# canvases.append(c2)
# ROOT.gROOT.GetListOfCanvases().Draw()
Position of 244456 on the distirbution¶
In [3]:
def read_branch(datafiles, jet_type, var, entrystop=None, entrystart=None, verbose=False, apply_func=None, cache=None):
if hasattr(datafiles, '__iter__') and type(datafiles) != str:
# multiple datafiles
vals = np.array([])
for f in datafiles:
if verbose: print(f)
kwargs = dict(datafiles=f, jet_type=jet_type, var=var, entrystop=entrystop, entrystart=entrystart, verbose=verbose, apply_func=apply_func, cache=cache)
v = read_branch(**kwargs)
vals = np.hstack([vals, v])
return vals
# single datafile
froot = uproot.open(datafiles)
if any(['JetPY' in str(k) for k in froot.keys()]):
tree_name_core = 'JetTree_AliAnalysisTaskJetExtractor_JetPY_AKTChargedR040_tracks_pT0150_E_scheme_'
else:
tree_name_core = 'JetTree_AliAnalysisTaskJetExtractor_Jet_AKTChargedR040_tracks_pT0150_E_scheme_'
if apply_func:
v0 = froot[tree_name_core+jet_type].array(var, flatten=False, entrystop=entrystop, entrystart=entrystart, cache=cache)
vals = [apply_func(arr) if len(arr) else 0 for arr in v0]
else:
vals = froot[tree_name_core+jet_type].array(var, flatten=True , entrystop=entrystop, entrystart=entrystart, cache=cache)
return vals
In [9]:
data_files = glob.glob('../../HF-jets/ana_results/iter2/LHC15n/myOutputDir/*/AnalysisResults.root')
for var in ['Event_Vertex_X', 'Event_Vertex_Y', 'Event_Multiplicity',
# 'Event_BackgroundDensity', 'Event_BackgroundDensityMass',
'Jet_Pt',
'Jet_Phi', 'Jet_Eta',
'Jet_Area',
'Jet_NumTracks',
'Jet_NumSecVertices',]:
print(var)
means, stds = [], []
for f in data_files:
vals = read_branch(f, 'allJets', var)
mean, std = np.mean(vals), np.std(vals)
means.append(mean)
stds.append(std)
# print(f'{f.split("000")[1][:6]}: {mean:.4f}, {std:.4f}')
if '244456' in f: mean244456, std244456 = mean, std
if '244453' in f: mean244453, std244453 = mean, std
fig,axes = plt.subplots(ncols=2, figsize=(10,5))
axes[0].hist(means, bins=15, histtype='step', linewidth=2)
axes[1].hist(stds, bins=15, histtype='step', linewidth=2)
axes[0].set_title('mean')
axes[1].set_title('std')
fig.suptitle(var)
xlim, ylim = axes[0].get_xlim(), axes[0].get_ylim()
deltax, deltay = xlim[1] - xlim[0] , ylim[1] - ylim[0]
# axes[0].set_xlim(xlim[0]-0.2*deltax, xlim[1]+0.2*deltax)
axes[0].arrow(mean244456, ylim[1], 0, -0.4*deltay, width=0.02*deltax, head_width=0.05*deltax, head_length=0.05*deltay, color='red')
axes[0].arrow(mean244453, ylim[1], 0, -0.4*deltay, width=0.015*deltax, head_width=0.04*deltax, head_length=0.05*deltay, color='orange')
xlim, ylim = axes[1].get_xlim(), axes[1].get_ylim()
deltax, deltay = xlim[1] - xlim[0] , ylim[1] - ylim[0]
# axes[1].set_xlim(xlim[0]-0.2*deltax, xlim[1]+0.2*deltax)
axes[1].arrow(std244456, ylim[1], 0, -0.4*deltay, width=0.02*deltax, head_width=0.05*deltax, head_length=0.05*deltay, color='red')
axes[1].arrow(std244453, ylim[1], 0, -0.4*deltay, width=0.015*deltax, head_width=0.04*deltax, head_length=0.05*deltay, color='orange')
In [ ]: