Belle II Software  release-05-01-25
v0ValidationCreatePlots.py
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 
4 """
5 <header>
6  <contact>software-tracking@belle2.org</contact>
7  <input>V0ValidationHarvested.root</input>
8  <description>This module creates efficiency plots for the V0 validation.</description>
9 </header>
10 """
11 
12 
13 import numpy
14 import ROOT
15 
16 
18 
19  """Reads the output created by the V0Harvester and creates plots from it."""
20 
21  def __init__(self, input_file='../V0ValidationHarvested.root', output_file='V0Validation.root'):
22  """Reads the output created by the V0Harvester defines histograms which will be filled later.
23 
24  :param input_file: Output of V0ValidationHarvester.
25  :param output_file: Plots displayed in the V0Validation.
26  """
27 
28  self.input_file = input_file
29 
30  self.output_file = output_file
31 
32 
33  self.hist_r = ROOT.TH1F("", "True R", 20, 0, 20)
34 
35  self.hist_theta = ROOT.TH1F("", "True Theta", 26, 20, 150)
36 
37  self.hist_phi = ROOT.TH1F("", "True Phi", 36, -180, 180)
38 
39  self.hist_p = ROOT.TH1F("", "True P", 25, 0.0, 1.0)
40 
41 
42  self.hist_r_found = ROOT.TH1F("", "Found R", 20, 0, 20)
43 
44  self.hist_theta_found = ROOT.TH1F("", "Found Theta", 26, 20, 150)
45 
46  self.hist_phi_found = ROOT.TH1F("", "Found Phi", 36, -180, 180)
47 
48  self.hist_p_found = ROOT.TH1F("", "Found P", 25, 0.0, 1.0)
49 
50 
51  self.hist_invariant_mass = ROOT.TH1F("", "", 60, 0.470, 0.530)
52 
53  self.hist_invariant_mass_res = ROOT.TH1F("", "", 40, -0.02, 0.02)
54 
55 
56  self.hist_chi2 = ROOT.TH1F("", "", 50, 0, 50)
57 
58  self.hist_chi2_inside = ROOT.TH1F("", "", 50, 0, 50)
59 
60  self.hist_chi2_outside = ROOT.TH1F("", "", 50, 0, 50)
61 
62 
63  self.hist_mass_vs_mc_mass = ROOT.TH2F("", "", 80, 0, 0.8, 80, 0, 0.8)
64 
65 
66  self.hist_invariant_mass_residuum = ROOT.TH1F("", "", 60, -0.05, 0.05)
67 
68  self.hist_r_residuum = ROOT.TH1F("", "", 60, -0.1, 0.1)
69 
70  self.hist_theta_residuum = ROOT.TH1F("", "", 60, -0.1, 0.1)
71 
72  self.hist_phi_residuum = ROOT.TH1F("", "", 60, -0.1, 0.1)
73 
74  self.hist_p_residuum = ROOT.TH1F("", "", 60, -0.05, 0.05)
75 
76  def collect_histograms(self):
77  """Fills the defined histograms with the V0Harvester data.
78 
79  :return: self
80  """
81  input_root_file = ROOT.TFile.Open(self.input_file, "READ")
82 
83  for event in input_root_file.V0Harvester_tree:
84  self.hist_r.Fill(event.R_MC)
85  self.hist_theta.Fill(numpy.rad2deg(event.THETA_MC))
86  self.hist_phi.Fill(numpy.rad2deg(event.PHI_MC))
87  self.hist_p.Fill(event.P_MC)
88 
89  if event.FOUND:
90  self.hist_r_found.Fill(event.R_MC)
91  self.hist_theta_found.Fill(numpy.rad2deg(event.THETA_MC))
92  self.hist_phi_found.Fill(numpy.rad2deg(event.PHI_MC))
93  self.hist_p_found.Fill(event.P_MC)
94 
95  self.hist_invariant_mass.Fill(event.M)
96 
97  # Residuum histograms
98  self.hist_invariant_mass_residuum.Fill(event.M - event.M_MC)
99  self.hist_r_residuum.Fill(event.R - event.R_MC)
100  self.hist_theta_residuum.Fill(event.THETA - event.THETA_MC)
101  self.hist_phi_residuum.Fill(event.PHI - event.PHI_MC)
102  self.hist_p_residuum.Fill(event.P - event.P_MC)
103 
104  self.hist_chi2.Fill(event.CHI2)
105  if event.R_MC > 1.0:
106  self.hist_chi2_outside.Fill(event.CHI2)
107  else:
108  assert event.R_MC <= 1.0
109  self.hist_chi2_inside.Fill(event.CHI2)
110 
111  self.hist_mass_vs_mc_mass.Fill(event.M, event.M_MC)
112 
113  input_root_file.Close()
114  return self
115 
116  @staticmethod
117  def efficiency_plot(found, all, title, x_variable, x_unit, description='', check='', contact='', meta_options=''):
118  """Create an efficiency plot.
119 
120  :param found: Histogram with all found entries (i.e. reconstructed).
121  :param all: Histogram with all entries (i.e. MCTruth).
122  :param title: Title of the histogram.
123  :param x_variable: x variable.
124  :param x_unit: x unit.
125  :param description: Description text shown on the validation page.
126  :param check: Check text shown on the validation page.
127  :param contact: Contact text shown on the validation page.
128  :param meta_options: Meta options for the validation page.
129  :return: ROOT.TEfficiency
130  """
131  efficiency = ROOT.TEfficiency(found, all)
132  efficiency.SetName("".join(title.split()))
133  ylabel = 'Efficiency / ({} {})'.format((found.GetXaxis().GetXmax() -
134  found.GetXaxis().GetXmin()) / found.GetNbinsX(), x_unit)
135  efficiency.SetTitle("{};{} / ({});{}".format(title, x_variable, x_unit, ylabel))
136  efficiency.GetListOfFunctions().Add(ROOT.TNamed('Description', description))
137  efficiency.GetListOfFunctions().Add(ROOT.TNamed('Check', check))
138  efficiency.GetListOfFunctions().Add(ROOT.TNamed('Contact', contact))
139  efficiency.GetListOfFunctions().Add(ROOT.TNamed('MetaOptions', meta_options))
140  return efficiency
141 
142  @staticmethod
143  def histogram_plot(hist, title, x_variable, x_unit=None, description='', check='', contact='', meta_options=''):
144  """Create (annotate) an histogram plot.
145 
146  :param hist: TH1F
147  :param title: Title of the histogram.
148  :param x_variable: x variable.
149  :param x_unit: x unit.
150  :param description: Description text shown on the validation page.
151  :param check: Check text shown on the validation page.
152  :param contact: Contact text shown on the validation page.
153  :param meta_options: Meta options for the validation page.
154  :return: modified hist
155  """
156  hist.SetName("".join(title.split()))
157  xlabel = '{} / ({})'.format(x_variable, x_unit) if x_unit is not None else '{}'.format(x_variable)
158  ylabel = 'Entries / ({} {})'.format((hist.GetXaxis().GetXmax() -
159  hist.GetXaxis().GetXmin()) /
160  hist.GetNbinsX(), x_unit) if x_unit is not None \
161  else 'Entries / ({})'.format((hist.GetXaxis().GetXmax() - hist.GetXaxis().GetXmin()) / hist.GetNbinsX())
162  hist.SetTitle("{};{};{}".format(title, xlabel, ylabel))
163  hist.GetListOfFunctions().Add(ROOT.TNamed('Description', description))
164  hist.GetListOfFunctions().Add(ROOT.TNamed('Check', check))
165  hist.GetListOfFunctions().Add(ROOT.TNamed('Contact', contact))
166  hist.GetListOfFunctions().Add(ROOT.TNamed('MetaOptions', meta_options))
167  return hist
168 
169  @staticmethod
170  def histogram_2d_plot(hist, title, x_variable, y_variable, x_unit=None, y_unit=None,
171  description='', check='', contact='', meta_options=''):
172  """Create a 2d hisogram plot.
173 
174  :param hist: TH2F
175  :param title: Title of the histogram.
176  :param x_variable: x variable
177  :param y_variable: y variable
178  :param x_unit: x unit
179  :param y_unit: y unit
180  :param description: Description text shown on the validation page.
181  :param check: Check text shown on the validation page.
182  :param contact: Contact text shown on the validation page.
183  :param meta_options: Meta options for the validation page.
184  :return: ROOT.TEfficiency
185  :return:
186  """
187  hist.SetName("".join(title.split()))
188  xlabel = '{} / ({})'.format(x_variable, x_unit) if x_unit is not None else '{}'.format(x_variable)
189  ylabel = '{} / ({})'.format(y_variable, y_unit) if y_unit is not None else '{}'.format(y_variable)
190  hist.SetTitle("{};{};{}".format(title, xlabel, ylabel))
191  hist.GetListOfFunctions().Add(ROOT.TNamed('Description', description))
192  hist.GetListOfFunctions().Add(ROOT.TNamed('Check', check))
193  hist.GetListOfFunctions().Add(ROOT.TNamed('Contact', contact))
194  hist.GetListOfFunctions().Add(ROOT.TNamed('MetaOptions', meta_options))
195  return hist
196 
197  def plot(self):
198  """Create plots with the data filled with 'collect_histograms'.
199 
200  :return: self
201  """
202  output_root_file = ROOT.TFile.Open(self.output_file, "RECREATE")
203 
204  V0ValidationPlots.efficiency_plot(self.hist_r_found, self.hist_r, 'Efficiency vs R', 'r', 'cm',
205  description='Reconstruction Efficiency vs. r (perpendicular)',
206  check='',
207  contact='software-tracking@belle2.org',
208  meta_options='').Write()
209 
210  V0ValidationPlots.efficiency_plot(self.hist_theta_found, self.hist_theta, 'Efficiency vs Theta', 'Theta', 'deg',
211  description='Reconstruction Efficiency vs. theta',
212  check='',
213  contact='software-tracking@belle2.org',
214  meta_options='').Write()
215 
216  V0ValidationPlots.efficiency_plot(self.hist_phi_found, self.hist_phi, 'Efficiency vs Phi', 'Phi', 'deg',
217  description='Reconstruction Efficiency vs phi',
218  check='',
219  contact='software-tracking@belle2.org',
220  meta_options='').Write()
221 
222  V0ValidationPlots.efficiency_plot(self.hist_p_found, self.hist_p, 'Efficiency vs P', 'P', 'GeV',
223  description='Reconstruction Efficiency vs momentum',
224  check='',
225  contact='software-tracking@belle2.org',
226  meta_options='').Write()
227 
228  V0ValidationPlots.histogram_plot(self.hist_invariant_mass, "KShort Invariant Mass", "m", "GeV",
229  description='Reconstructed invariant mass of KShorts with standard reconstruction',
230  check='Invariant mass peak around KShort nominal mass 497.61 MeV.',
231  contact='software-tracking@belle2.org',
232  meta_options='').Write()
233 
234  V0ValidationPlots.histogram_plot(self.hist_invariant_mass_residuum, "KShort Invariant Mass Residuum", "Rec - MC", "GeV",
235  description='Invariant mass residuum of KShorts with standard reconstruction',
236  check='',
237  contact='software-tracking@belle2.org',
238  meta_options='').Write()
239 
240  V0ValidationPlots.histogram_plot(self.hist_chi2, "Chi2 of Vertex Fits.", "Chi2", None,
241  description='Chi2 distributions of the vertex fits.',
242  check='Check if distribution looks like a Chi2 distribution with 1 dof',
243  contact='software-tracking@belle2.org',
244  meta_options='').Write()
245 
246  V0ValidationPlots.histogram_plot(self.hist_chi2_inside, "Chi2 of Vertex Fits Inside Beampipe.", "Chi2", None,
247  description='Chi2 distributions of the vertex fits inside the beampipe.',
248  check='Check if distribution looks like a Chi2 distribution with 1 dof',
249  contact='software-tracking@belle2.org',
250  meta_options='expert').Write()
251 
252  V0ValidationPlots.histogram_plot(self.hist_chi2_outside, "Chi2 of Vertex Fits Outside Beampipe.", "Chi2", None,
253  description='Chi2 distributions of the vertex fits outside the beampipe.',
254  check='Check if distribution looks like a Chi2 distribution with 1 dof',
255  contact='software-tracking@belle2.org',
256  meta_options='expert').Write()
257 
258  V0ValidationPlots.histogram_2d_plot(self.hist_mass_vs_mc_mass, "Reconstructed vs MC Mass.",
259  "Reconstructed Mass", "GeV", "MC Mass", "GeV",
260  description="Reconstructed mass vs invariant Mass.",
261  check="",
262  contact="software-tracking@belle2.org",
263  meta_options='expert').Write()
264 
265  V0ValidationPlots.histogram_plot(self.hist_r_residuum, "KShort R Residuum", "Rec - MC", "cm",
266  description='R residuum of KShorts with standard reconstruction',
267  check='',
268  contact='software-tracking@belle2.org',
269  meta_options='').Write()
270  V0ValidationPlots.histogram_plot(self.hist_theta_residuum, "KShort Theta Residuum", "Rec - MC", "rad",
271  description='Theta residuum of KShorts with standard reconstruction',
272  check='',
273  contact='software-tracking@belle2.org',
274  meta_options='').Write()
275  V0ValidationPlots.histogram_plot(self.hist_phi_residuum, "KShort Phi Residuum", "Rec - MC", "rad",
276  description='Phi residuum of KShorts with standard reconstruction',
277  check='',
278  contact='software-tracking@belle2.org',
279  meta_options='').Write()
280  V0ValidationPlots.histogram_plot(self.hist_p_residuum, "KShort Momentum Residuum", "Rec - MC", "GeV",
281  description='Momentum residuum of KShorts with standard reconstruction',
282  check='',
283  contact='software-tracking@belle2.org',
284  meta_options='').Write()
285 
286  output_root_file.Write()
287  output_root_file.Close()
288  return self
289 
290 
291 if __name__ == '__main__':
v0ValidationCreatePlots.V0ValidationPlots.plot
def plot(self)
Definition: v0ValidationCreatePlots.py:197
v0ValidationCreatePlots.V0ValidationPlots.hist_invariant_mass_res
hist_invariant_mass_res
Invariant mass residual histogram.
Definition: v0ValidationCreatePlots.py:53
v0ValidationCreatePlots.V0ValidationPlots.output_file
output_file
Output file.
Definition: v0ValidationCreatePlots.py:30
v0ValidationCreatePlots.V0ValidationPlots.hist_theta_found
hist_theta_found
Reconstructed Theta.
Definition: v0ValidationCreatePlots.py:44
v0ValidationCreatePlots.V0ValidationPlots.histogram_2d_plot
def histogram_2d_plot(hist, title, x_variable, y_variable, x_unit=None, y_unit=None, description='', check='', contact='', meta_options='')
Definition: v0ValidationCreatePlots.py:170
plot
Definition: plot.py:1
v0ValidationCreatePlots.V0ValidationPlots.hist_chi2
hist_chi2
Chi2 of all vertex fits.
Definition: v0ValidationCreatePlots.py:56
v0ValidationCreatePlots.V0ValidationPlots
Definition: v0ValidationCreatePlots.py:17
v0ValidationCreatePlots.V0ValidationPlots.hist_chi2_outside
hist_chi2_outside
Chi2 of vertex fit outside beampipe.
Definition: v0ValidationCreatePlots.py:60
v0ValidationCreatePlots.V0ValidationPlots.hist_phi
hist_phi
MCTruth Phi.
Definition: v0ValidationCreatePlots.py:37
v0ValidationCreatePlots.V0ValidationPlots.collect_histograms
def collect_histograms(self)
Definition: v0ValidationCreatePlots.py:76
v0ValidationCreatePlots.V0ValidationPlots.hist_p_residuum
hist_p_residuum
Momentum residuum histogram.
Definition: v0ValidationCreatePlots.py:74
v0ValidationCreatePlots.V0ValidationPlots.hist_phi_residuum
hist_phi_residuum
Phi residuum histogram.
Definition: v0ValidationCreatePlots.py:72
v0ValidationCreatePlots.V0ValidationPlots.hist_p_found
hist_p_found
Reconstructed Momentum.
Definition: v0ValidationCreatePlots.py:48
v0ValidationCreatePlots.V0ValidationPlots.hist_chi2_inside
hist_chi2_inside
Chi2 of vertex fit inside beampipe.
Definition: v0ValidationCreatePlots.py:58
v0ValidationCreatePlots.V0ValidationPlots.hist_r_found
hist_r_found
Reconstructed R.
Definition: v0ValidationCreatePlots.py:42
v0ValidationCreatePlots.V0ValidationPlots.hist_invariant_mass
hist_invariant_mass
Invariant mass histogram.
Definition: v0ValidationCreatePlots.py:51
v0ValidationCreatePlots.V0ValidationPlots.hist_r
hist_r
MCTruth R.
Definition: v0ValidationCreatePlots.py:33
v0ValidationCreatePlots.V0ValidationPlots.hist_phi_found
hist_phi_found
Reconstructed Phi.
Definition: v0ValidationCreatePlots.py:46
v0ValidationCreatePlots.V0ValidationPlots.hist_theta
hist_theta
MCTruth Theta.
Definition: v0ValidationCreatePlots.py:35
v0ValidationCreatePlots.V0ValidationPlots.efficiency_plot
def efficiency_plot(found, all, title, x_variable, x_unit, description='', check='', contact='', meta_options='')
Definition: v0ValidationCreatePlots.py:117
v0ValidationCreatePlots.V0ValidationPlots.hist_p
hist_p
MCTruth Momentum.
Definition: v0ValidationCreatePlots.py:39
v0ValidationCreatePlots.V0ValidationPlots.histogram_plot
def histogram_plot(hist, title, x_variable, x_unit=None, description='', check='', contact='', meta_options='')
Definition: v0ValidationCreatePlots.py:143
v0ValidationCreatePlots.V0ValidationPlots.__init__
def __init__(self, input_file='../V0ValidationHarvested.root', output_file='V0Validation.root')
Definition: v0ValidationCreatePlots.py:21
v0ValidationCreatePlots.V0ValidationPlots.input_file
input_file
Input file.
Definition: v0ValidationCreatePlots.py:28
v0ValidationCreatePlots.V0ValidationPlots.hist_theta_residuum
hist_theta_residuum
Theta residuum histogram.
Definition: v0ValidationCreatePlots.py:70
v0ValidationCreatePlots.V0ValidationPlots.hist_r_residuum
hist_r_residuum
R residuum histogram.
Definition: v0ValidationCreatePlots.py:68
v0ValidationCreatePlots.V0ValidationPlots.hist_mass_vs_mc_mass
hist_mass_vs_mc_mass
2D histogram; invariant mass vs reconstructed mass
Definition: v0ValidationCreatePlots.py:63
v0ValidationCreatePlots.V0ValidationPlots.hist_invariant_mass_residuum
hist_invariant_mass_residuum
Invariant mass residuum histogram.
Definition: v0ValidationCreatePlots.py:66