Belle II Software  release-05-01-25
inputVariablesPlots.py
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 
4 
17 
18 import os
19 import glob
20 import sys
21 import math
22 from matplotlib.ticker import FormatStrFormatter
23 import matplotlib.pyplot as plt
24 import ROOT
25 from ROOT import Belle2
26 import flavorTagger as ft
27 import basf2_mva
28 from defaultEvaluationParameters import r_subsample, r_size, categories
29 from array import array
30 
31 import numpy as np
32 import matplotlib as mpl
33 mpl.use('Agg')
34 mpl.rcParams.update({'font.size': 22})
35 mpl.rcParams['text.usetex'] = True
36 mpl.rcParams['text.latex.preamble'] = [r"\usepackage{amsmath}"]
37 
38 ROOT.gROOT.SetBatch(True)
39 
40 if len(sys.argv) != 5:
41  sys.exit("Must provide 4 arguments: [Belle or Belle2] [BGx0 or BGx1] [training] [workingDirectory]"
42  )
43 
44 belleOrBelle2 = sys.argv[1]
45 MCtype = str(sys.argv[2])
46 decayChannelTrainedOn = str(sys.argv[3]) # Decay channel of the weight files "JPsiKs" or "nunubar"
47 workingDirectory = sys.argv[4]
48 
49 filesDirectory = workingDirectory + '/FlavorTagging/TrainedMethods'
50 
51 if decayChannelTrainedOn == 'JPsiKs':
52  decayChannelTrainedOn = 'JpsiKs_mu'
53 
54 weightFiles = 'B2' + decayChannelTrainedOn + MCtype
55 
56 
57 ROOT.TH1.SetDefaultSumw2()
58 
59 allInputVariables = []
60 
61 
62 ft.WhichCategories(categories)
63 ft.setVariables()
64 
65 belleOrBelle2Flag = belleOrBelle2
66 
67 identifiersExtraInfosDict = dict()
68 identifiersExtraInfosKaonPion = []
69 
70 if belleOrBelle2 == "Belle":
71  unitImp = "cm"
72 
73 dBw = 50
74 
75 pBins = 50
76 fBins = 100
77 
78 unitImp = "mm"
79 
80 variablesPlotParamsDict = {
81  'useCMSFrame(p)': ['useCMSFrame__bop__bc', pBins, 0, 3, r'$p^*\ [{\rm GeV}/c]$', r"{\rm GeV}/c\, "],
82  'useCMSFrame(pt)': ['useCMSFrame__bopt__bc', pBins, 0, 3, r'$p_{\rm t}^*\ [{\rm GeV}/c]$', r"{\rm GeV}/c\, "],
83  'p': ['p', pBins, 0, 3, r'$p\ [{\rm GeV}/c]$', r"{\rm GeV}/c\, "],
84  'pt': ['pt', pBins, 0, 3, r'$p_{\rm t}\ [{\rm GeV}/c]$', r"{\rm GeV}/c\, "],
85  'cosTheta': ['cosTheta', dBw, -1, 1.01, r'$\cos{\theta}$', ""],
86  ft.eId[ft.getBelleOrBelle2()]: [Belle2.makeROOTCompatible(ft.eId[ft.getBelleOrBelle2()]),
87  dBw, 0, 1.01, r'$\mathcal{L}_{e}$', ""],
88  'eid_dEdx': ['eid_dEdx', dBw, 0, 1.01, r'$\mathcal{L}_{e}^{{\rm d}E/{\rm d}x}$', ""],
89  'eid_TOP': ['eid_TOP', dBw, 0, 1.01, r'$\mathcal{L}_{e}^{\rm TOP}$', ""],
90  'eid_ARICH': ['eid_ARICH', dBw, 0, 1.01, r'$\mathcal{L}_{e}^{\rm ARICH}$', ""],
91  'eid_ECL': ['eid_ECL', dBw, 0, 1.01, r'$\mathcal{L}_{e}^{\rm ECL}$', ""],
92  'BtagToWBosonVariables(recoilMassSqrd)': ['BtagToWBosonVariables__borecoilMassSqrd__bc', fBins,
93  0, 100, r'$M_{\rm rec}^2\ [{\rm GeV}^2/c^4]$', r"{\rm GeV}^2/c^4"],
94  'BtagToWBosonVariables(pMissCMS)': [
95  'BtagToWBosonVariables__bopMissCMS__bc', 60, 0, 3.6,
96  r'$p^*_{\rm miss}\ [{\rm GeV}/c]$', r"{\rm GeV}/c\, "],
97  'BtagToWBosonVariables(cosThetaMissCMS)': ['BtagToWBosonVariables__bocosThetaMissCMS__bc',
98  dBw, -1, 1.01, r'$\cos{\theta^*_{\rm miss}}$', ""],
99  'BtagToWBosonVariables(EW90)': ['BtagToWBosonVariables__boEW90__bc',
100  dBw, 0, 4, r'$E_{90}^{W}\ [{\rm GeV}]$', r"{\rm GeV}\, "],
101  'BtagToWBosonVariables(recoilMass)': ['sqrt(abs(BtagToWBosonVariables__borecoilMassSqrd__bc))',
102  dBw, 0, 12, r'$M_{\rm rec}\ [{\rm GeV}/c^2]$', r"{\rm GeV}/c^2\, "],
103  'cosTPTO': ['cosTPTO', dBw, 0, 1.01, r'$\vert\cos{\theta^*_{\rm T}}\vert$', ""],
104  'ImpactXY': ['ImpactXY', dBw, 0, 0.5, r'$d_0\ [{\rm ' + unitImp + '}]$', r"{\rm " + unitImp + r"}\, "],
105  'z0': ['z0', dBw, 0, 1.0, r'$z_0\ [{\rm ' + unitImp + '}]$', r"{\rm " + unitImp + r"}\, "],
106  'y': ['y', dBw, -0.15, 0.15, r'$y_0\ [{\rm ' + unitImp + '}]$', r"{\rm " + unitImp + r"}\, "],
107  'OBoost': ['OBoost', dBw, -0.15, 0.15, r'$d_0^\prime\ [{\rm ' + unitImp + '}]$', r"{\rm " + unitImp + r"}\, "],
108  'distance': ['distance', dBw, 0, 1.5, r'$\xi_0\ [{\rm ' + unitImp + '}]$', r"{\rm " + unitImp + r"}\, "],
109  'chiProb': ['chiProb', dBw, 0, 1.01, r'$p$-${\rm value}$', ""],
110  ft.muId[ft.getBelleOrBelle2()]: [Belle2.makeROOTCompatible(ft.muId[ft.getBelleOrBelle2()]),
111  dBw, 0, 1.01, r'$\mathcal{L}_{\mu}$', ""],
112  'muid_dEdx': ['muid_dEdx', dBw, 0, 1.01, r'$\mathcal{L}_{\mu}^{{\rm d}E/{\rm d}x}$', ""],
113  'muid_TOP': ['muid_TOP', dBw, 0, 1.01, r'$\mathcal{L}_{\mu}^{\rm TOP}$', ""],
114  'muid_ARICH': ['muid_ARICH', dBw, 0, 1.01, r'$\mathcal{L}_{\mu}^{\rm ARICH}$', ""],
115  'muid_KLM': ['muid_KLM', dBw, 0, 1.01, r'$\mathcal{L}_{\mu}^{\rm KLM}$', ""],
116  ft.KId[ft.getBelleOrBelle2()]: [Belle2.makeROOTCompatible(ft.KId[ft.getBelleOrBelle2()]),
117  dBw, 0, 1.01, r'$\mathcal{L}_{K}$', ""],
118  'Kid_dEdx': ['Kid_dEdx', dBw, 0, 1.01, r'$\mathcal{L}_{K}^{{\rm d}E/{\rm d}x}$', ""],
119  'Kid_TOP': ['Kid_TOP', dBw, 0, 1.01, r'$\mathcal{L}_{K}^{\rm TOP}$', ""],
120  'Kid_ARICH': ['Kid_ARICH', dBw, 0, 1.01, r'$\mathcal{L}_{K}^{\rm ARICH}$', ""],
121  'NumberOfKShortsInRoe': ['NumberOfKShortsInRoe', dBw, 0, 12, r'$n_{K^0_S}$', ""],
122  'ptTracksRoe': ['ptTracksRoe', fBins, 0, 6, r'$\Sigma\, p_{\rm t}^2\ [{\rm GeV^2}/c^2]$',
123  r"{\rm GeV^2}/c^2"],
124  'extraInfo(isRightCategory(Kaon))': ['extraInfo__boisRightCategory__boKaon__bc__bc',
125  dBw, 0, 1.01, r"$y_{\rm Kaon}$", ""],
126  'HighestProbInCat(pi+:inRoe, isRightCategory(SlowPion))': [
127  'HighestProbInCat__bopi__pl__clinRoe__cm__spisRightCategory__boSlowPion__bc__bc',
128  dBw, 0, 1.01, r"$y_{\rm SlowPion}$", ""],
129  'KaonPionVariables(cosKaonPion)': ['KaonPionVariables__bocosKaonPion__bc',
130  dBw, -1, 1.01, r'$\cos{\theta^*_{K\pi}}$', ""],
131  'KaonPionVariables(HaveOpositeCharges)': ['KaonPionVariables__boHaveOpositeCharges__bc',
132  dBw, 0, 1.01, r'$\frac{1 - q_{K} \cdot q_\pi}{2}$', ""],
133  'pionID': ['pionID', dBw, 0, 1.01, r'$\mathcal{L}_{\pi}$', ""],
134  'piid_dEdx': ['piid_dEdx', dBw, 0, 1.01, r'$\mathcal{L}_{\pi}^{{\rm d}E/{\rm d}x}$', ""],
135  'piid_TOP': ['piid_TOP', dBw, 0, 1.01, r'$\mathcal{L}_{\pi}^{\rm TOP}$', ""],
136  'piid_ARICH': ['piid_ARICH', dBw, 0, 1.01, r'$\mathcal{L}_{\pi}^{\rm ARICH}$', ""],
137  'pi_vs_edEdxid': ['pi_vs_edEdxid', dBw, 0, 1.01, r'$\mathcal{L}_{\pi/e}^{{\rm d}E/{\rm d}x}$', ""],
138  'FSCVariables(pFastCMS)': ['FSCVariables__bopFastCMS__bc', pBins, 0, 3, r'$p^*_{\rm Fast}\ [{\rm GeV}/c]$', r"{\rm GeV}/c\, "],
139  'FSCVariables(cosSlowFast)': ['FSCVariables__bocosSlowFast__bc', dBw, -1, 1.01, r'$\cos{\theta^*_{\rm SlowFast}}$', ''],
140  'FSCVariables(cosTPTOFast)': ['FSCVariables__bocosTPTOFast__bc', dBw, 0, 1.01, r'$\vert\cos{\theta^*_{\rm T, Fast}}\vert$', ''],
141  'FSCVariables(SlowFastHaveOpositeCharges)': ['FSCVariables__boSlowFastHaveOpositeCharges__bc',
142  dBw, 0, 1.01, r'$\frac{1 - q_{\rm Slow} \cdot q_{\rm Fast}}{2}$', ""],
143  'lambdaFlavor': ['lambdaFlavor', dBw, -1, 1.01, r'$q_{\Lambda}$', ""],
144  'M': ['M', dBw, 1.08, 1.22, r'$M_{\Lambda}\ [{\rm GeV}/c^2]$', r"{\rm MeV}/c^2\, "],
145  'cosAngleBetweenMomentumAndVertexVector': ['cosAngleBetweenMomentumAndVertexVector',
146  dBw, -1, 1.01,
147  r'$\cos{\theta_{\boldsymbol{x}_{\Lambda},\boldsymbol{p}_{\Lambda}}}$', ""],
148  'lambdaZError': ['lambdaZError', dBw, 0, 0.05, r'$\sigma_{\Lambda}^{zz}$', r"{\rm mm}\, "],
149  'daughter(0,p)': ['daughter__bo0__cmp__bc', dBw, 0, 1, r'$p_{\pi}\ [{\rm GeV}/c]$', r"{\rm GeV}/c\, "],
150  'daughter(0,useCMSFrame(p))': ['daughter__bo0__cmuseCMSFrame__bop__bc__bc',
151  dBw, 0, 1, r'$p^*_{\pi}\ [{\rm GeV}/c]$', r"{\rm GeV}/c\, "],
152  'daughter(1,p)': ['daughter__bo1__cmp__bc', dBw, 0, 2, r'$p_{p}\ [{\rm GeV}/c]$', r"{\rm GeV}/c"],
153  'daughter(1,useCMSFrame(p))': ['daughter__bo1__cmuseCMSFrame__bop__bc__bc',
154  dBw, 0, 2, r'$p^*_{p}\ [{\rm GeV}/c]$', r"{\rm GeV}/c\, "],
155  'daughter(1,protonID)': ['daughter__bo1__cmprotonID__bc', dBw, 0, 1.01, r'$\mathcal{L}_{p}$', ""],
156  'daughter(0,pionID)': ['daughter__bo0__cmpionID__bc', dBw, 0, 1.01, r'$\mathcal{L}_{\pi}$', ""]}
157 
158 if not Belle2.FileSystem.findFile('./InputVariablesPlots', True):
159  os.mkdir('./InputVariablesPlots')
160 
161 
162 def plotInputVariablesOfFlavorTagger():
163  """
164  Makes plots of the distribution of the input variables of the flavor tagger
165  for each category distinguishing between the target particles of the category (signal)
166  and all the other (bkg.)
167  """
168 
169  for (particleList, category, combinerVariable) in ft.eventLevelParticleLists:
170 
171  # if category != "SlowPion":
172  # continue
173 
174  if not Belle2.FileSystem.findFile('./InputVariablesPlots/' + category, True):
175  os.mkdir('./InputVariablesPlots/' + category)
176 
177  if particleList not in identifiersExtraInfosDict and category != 'KaonPion':
178  identifiersExtraInfosDict[particleList] = []
179 
180  methodPrefixEventLevel = "FlavorTagger_" + belleOrBelle2Flag + "_" + weightFiles + 'EventLevel' + category + 'FBDT'
181  treeName = methodPrefixEventLevel + "_tree"
182  targetVariable = 'isRightCategory(' + category + ')'
183 
184  tree = ROOT.TChain(treeName)
185 
186  workingFiles = glob.glob(filesDirectory + '/' + methodPrefixEventLevel + 'sampled*.root')
187  # print("workingFiles = ", workingFiles)
188  # workingFiles = glob.glob(filesDirectory + '/' + methodPrefixEventLevel + 'sampled1?.root')
189 
190  for iFile in workingFiles:
191  # if Belle2.FileSystem.findFile(workingFile):
192  tree.AddFile(iFile)
193 
194  categoryInputVariables = []
195  trulyUsedInputVariables = []
196  for iVariable in tree.GetListOfBranches():
197 
198  managerVariableName = str(Belle2.invertMakeROOTCompatible(iVariable.GetName()))
199 
200  if managerVariableName in ft.variables[category] or managerVariableName == 'distance' or \
201  managerVariableName == 'z0' or managerVariableName == 'ImpactXY' or \
202  managerVariableName == 'y' or managerVariableName == 'OBoost':
203  if managerVariableName in categoryInputVariables:
204  continue
205 
206  categoryInputVariables.append(managerVariableName)
207  if managerVariableName in ft.variables[category]:
208  allInputVariables.append((category, managerVariableName))
209  trulyUsedInputVariables.append((category, managerVariableName))
210 
211  if managerVariableName not in identifiersExtraInfosDict[particleList] and category != 'KaonPion':
212  identifiersExtraInfosDict[particleList].append(managerVariableName)
213 
214  elif category == 'KaonPion' and managerVariableName not in identifiersExtraInfosKaonPion:
215  identifiersExtraInfosKaonPion.append(managerVariableName)
216 
217  # if managerVariableName not in variablesPlotParamsDict:
218  # variablesPlotParamsDict[managerVariableName] =
219  # [iVariable.GetName(), 100, 0, 2, iVariable.GetName(), "unit"]
220 
221  print("The number of variables used in " + category + " is = ", len(trulyUsedInputVariables))
222 
223  if category != 'KaonPion' and category != 'FSC' and category != 'MaximumPstar' and \
224  category != 'FastHadron' and category != 'Lambda':
225  categoryInputVariables.append('BtagToWBosonVariables(recoilMass)')
226 
227  for inputVariable in categoryInputVariables:
228 
229  print(inputVariable)
230 
231  nBins = variablesPlotParamsDict[inputVariable][1]
232  limXmin = variablesPlotParamsDict[inputVariable][2]
233  limXmax = variablesPlotParamsDict[inputVariable][3]
234 
235  if category == "SlowPion":
236  if inputVariable == 'p' or inputVariable == 'useCMSFrame(p)' or \
237  inputVariable == 'pt' or inputVariable == 'useCMSFrame(pt)':
238  nBins = 150
239  limXmax = 1.5
240 
241  if inputVariable == 'distance':
242  nBins = 80
243  limXmax = 2.4
244 
245  if inputVariable == 'ImpactXY':
246  nBins = 80
247  limXmax = 0.8
248 
249  if category == "Lambda":
250  if inputVariable == 'distance':
251  # nBins = 25
252  limXmax = 10
253 
254  if inputVariable == 'chiProb':
255  nBins = 25
256 
257  signalHistogram = ROOT.TH1F("signal" + category + str(Belle2.makeROOTCompatible(inputVariable)), "",
258  nBins,
259  limXmin,
260  limXmax)
261  backgroundHistogram = ROOT.TH1F("bkg" + category + str(Belle2.makeROOTCompatible(inputVariable)), "",
262  nBins,
263  limXmin,
264  limXmax)
265 
266  factorMultiplication = str()
267 
268  if belleOrBelle2 == "Belle2" and ((category != "Lambda" and inputVariable == 'distance') or inputVariable ==
269  'z0' or inputVariable == 'ImpactXY' or inputVariable ==
270  'y' or inputVariable == 'OBoost'):
271  factorMultiplication = "*10 "
272 
273  tree.Draw(variablesPlotParamsDict[inputVariable][0] + factorMultiplication + ">> signal" + category +
274  str(Belle2.makeROOTCompatible(inputVariable)), Belle2.makeROOTCompatible(targetVariable) + " > 0")
275 
276  tree.Draw(variablesPlotParamsDict[inputVariable][0] + factorMultiplication + ">> bkg" + category +
277  str(Belle2.makeROOTCompatible(inputVariable)), Belle2.makeROOTCompatible(targetVariable) + " < 1")
278 
279  signalScalingFactor = signalHistogram.Integral()
280  backgroundScalingFactor = backgroundHistogram.Integral()
281 
282  if signalScalingFactor == 0:
283  signalScalingFactor = 1
284 
285  if backgroundScalingFactor == 0:
286  backgroundScalingFactor = 1
287 
288  signalHistogram.Scale(1 / signalScalingFactor)
289  backgroundHistogram.Scale(1 / backgroundScalingFactor)
290 
291  signalArray = np.zeros((signalHistogram.GetNbinsX(), 2))
292  backgroundArray = np.zeros((backgroundHistogram.GetNbinsX(), 2))
293 
294  for i in range(0, signalHistogram.GetNbinsX()):
295  signalArray[i] = np.array([signalHistogram.GetBinCenter(i + 1), signalHistogram.GetBinContent(i + 1)])
296  backgroundArray[i] = np.array([backgroundHistogram.GetBinCenter(i + 1), backgroundHistogram.GetBinContent(i + 1)])
297 
298  fig1 = plt.figure(1, figsize=(11, 10))
299 
300  # if inputVariable == 'Kid_dEdx' or inputVariable == 'muid_dEdx':
301  # ax1 = plt.axes([0.18, 0.17, 0.75, 0.8])
302  # if inputVariable == 'pi_vs_edEdxid':
303  # ax1 = plt.axes([0.18, 0.187, 0.75, 0.805])
304  # else:
305  ax1 = plt.axes([0.18, 0.2, 0.76, 0.705])
306 
307  # print(signalArray.shape, signalHistogram.GetNbinsX(), )
308  # print(signalArray)
309  ax1.hist(
310  signalArray[:, 0], weights=signalArray[:, 1], bins=signalHistogram.GetNbinsX(),
311  histtype='step',
312  edgecolor='r',
313  linewidth=4,
314  alpha=0.9,
315  label=r'${\rm Signal}$')
316 
317  ax1.hist(backgroundArray[:, 0], weights=backgroundArray[:, 1], bins=backgroundHistogram.GetNbinsX(),
318  histtype='step',
319  edgecolor='b', linewidth=4.5, linestyle='dashed', label=r'${\rm Background}$') # hatch='.',
320 
321  p1, = ax1.plot([], label=r'${\rm Signal}$', linewidth=5, linestyle='solid', alpha=0.9, c='r')
322  p2, = ax1.plot([], label=r'${\rm Background}$', linewidth=5.5, linestyle='dashed', c='b')
323 
324  binWidth = signalHistogram.GetBinWidth(2)
325 
326  if inputVariable == 'lambdaZError': # or inputVariable == 'ImpactXY' or\
327  # (category != "Lambda" and inputVariable == 'distance'):
328  binWidth = binWidth * 10
329 
330  if inputVariable == 'M':
331  binWidth = binWidth * 1000
332 
333  if category == "Lambda" and inputVariable == 'distance':
334  variablesPlotParamsDict[inputVariable][5] = r"{\rm cm}\, "
335 
336  binWidth = '{:8.2f}'.format(binWidth)
337 
338  xLabel = variablesPlotParamsDict[inputVariable][4]
339  legendLocation = 1
340 
341  if category == "FSC":
342  if inputVariable == 'cosTPTO':
343  xLabel = r'$\vert\cos{\theta^*_{\rm T, Slow}}\vert$'
344  if inputVariable == 'useCMSFrame(p)':
345  xLabel = r'$p^*_{\rm Slow}\ [{\rm GeV}/c]$'
346  if category == 'Lambda':
347  if inputVariable == 'useCMSFrame(p)':
348  xLabel = r'$p^*_{\Lambda}\ [{\rm GeV}/c]$'
349  if inputVariable == 'p':
350  xLabel = r'$p_{\Lambda}\ [{\rm GeV}/c]$'
351 
352  ax1.set_ylabel(r'${\rm Fraction\hspace{0.25em} of\hspace{0.25em} Events}\, /\, (\, ' + binWidth + r'\, ' +
353  variablesPlotParamsDict[inputVariable][5] + r')$', fontsize=46)
354  ax1.set_xlabel(xLabel, fontsize=65)
355  # plt.xticks([-1, -0.75, -0.5, -0.25, 0, 0.25, 0.5, 0.75, 1],
356  # [r'$-1$', r'', r'$-0.5$', r'', r'$0$', r'', r'$0.5$', r'', r'$1$'], rotation=0, size=40)
357  if inputVariable == 'extraInfo(isRightCategory(Kaon))' or \
358  inputVariable == 'HighestProbInCat(pi+:inRoe, isRightCategory(SlowPion))':
359  legendLocation = 3
360  ax1.set_yscale('log', nonposy='clip')
361  else:
362  ax1.yaxis.set_major_formatter(FormatStrFormatter(r'$%.2f$'))
363 
364  ax1.tick_params(axis='x', labelsize=50)
365  ax1.tick_params(axis='y', labelsize=40)
366 
367  if inputVariable == 'pi_vs_edEdxid':
368  ax1.xaxis.labelpad = 5
369  else:
370  ax1.xaxis.labelpad = 15
371 
372  if inputVariable.find('ARICH') != -1 or inputVariable.find('TOP') != -1 or \
373  inputVariable == 'cosTPTO' or inputVariable.find('KLM') != -1 or \
374  inputVariable == 'cosTheta' or inputVariable == 'FSCVariables(cosTPTOFast)' or \
375  inputVariable == 'KaonPionVariables(cosKaonPion)' or \
376  inputVariable == 'BtagToWBosonVariables(recoilMass)':
377  legendLocation = 2
378 
379  elif inputVariable == 'FSCVariables(SlowFastHaveOpositeCharges)' or \
380  inputVariable == 'KaonPionVariables(HaveOpositeCharges)' or inputVariable == "eid_ECL" or \
381  inputVariable.find('ID') != -1 or inputVariable.find('dEdx') != -1:
382  legendLocation = 9
383 
384  if inputVariable == 'muid_dEdx':
385  if category != 'KinLepton':
386  legendLocation = 8
387 
388  ax1.legend([p1, p2], [r'${\rm Signal}$', r'${\rm Bkgr.}$'], prop={
389  'size': 50}, loc=legendLocation, numpoints=1, handlelength=1)
390  ax1.grid(True)
391  # ax1.set_ylim(0, 1.4)
392  ax1.set_xlim(limXmin, limXmax)
393  plt.savefig('./InputVariablesPlots/' + category + '/' + category +
394  "_" + str(Belle2.makeROOTCompatible(inputVariable)) + '.pdf')
395  fig1.clear()
396 
397  signalHistogram.Delete()
398  backgroundHistogram.Delete()
399 
400 
401 if __name__ == '__main__':
402 
403  plotInputVariablesOfFlavorTagger()
404 
405  totalNumberOfVariables = 0
406 
407  for category in ft.variables:
408  totalNumberOfVariables += len(ft.variables[category])
409 
410  print("Total number of variables = ", totalNumberOfVariables)
411 
412  totalNumberOfCalculatedVariables = len(identifiersExtraInfosKaonPion)
413 
414  print("Calculations for Kaon-Pion Category = ", totalNumberOfCalculatedVariables)
415 
416  print("Variables per particle list:")
417  for particleList in identifiersExtraInfosDict:
418  print(particleList)
419  print(identifiersExtraInfosDict[particleList])
420  totalNumberOfCalculatedVariables += len(identifiersExtraInfosDict[particleList])
421 
422  print("Total number of calculated variables = ", totalNumberOfCalculatedVariables)
Belle2::invertMakeROOTCompatible
std::string invertMakeROOTCompatible(std::string str)
Invert makeROOTCompatible operation.
Definition: MakeROOTCompatible.cc:89
Belle2::makeROOTCompatible
std::string makeROOTCompatible(std::string str)
Remove special characters that ROOT dislikes in branch names, e.g.
Definition: MakeROOTCompatible.cc:74
Belle2::FileSystem::findFile
static std::string findFile(const std::string &path, bool silent=false)
Search for given file or directory in local or central release directory, and return absolute path if...
Definition: FileSystem.cc:147