Belle II Software  release-05-01-25
extmuid_createPlots.py
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 
4 
13 
14 """
15 <header>
16  <input>muon-ExtMuidValidation.root,
17  pion-ExtMuidValidation.root</input>
18  <contact>piilonen@vt.edu</contact>
19  <description>Create validation plots for Ext and Muid</description>
20 </header>
21 """
22 
23 import ROOT
24 
25 ROOT.PyConfig.IgnoreCommandLineOptions = True
26 ROOT.gROOT.SetBatch(True)
27 
28 from ROOT import TFile, TChain, TTree, TH1F, TH2F, TCanvas, TNamed, gStyle
29 import sys
30 import math
31 import numpy as np
32 from optparse import OptionParser
33 
34 # contact person information
35 # is added to the plot descriptions
36 CONTACT_PERSON = {'Name': 'Leo Piilonen',
37  'Email': 'piilonen@vt.edu'}
38 
39 
40 def main():
41  """Create validation plots for Ext and Muid"""
42 
43  option_parser = OptionParser()
44  option_parser.add_option('-i', '--input-file', dest='input_file',
45  default='../*on-ExtMuidValidation.root',
46  help='Root file with Ext and Muid validation data.'
47  )
48  option_parser.add_option('-o', '--output-file', dest='output_file',
49  default='ExtMuidValidation.root',
50  help='Root file with Ext and Muid validation histograms.')
51 
52  (options, args) = option_parser.parse_args()
53 
54  # load chain of input files
55  file_chain = TChain('tree')
56  file_chain.Add(options.input_file)
57 
58  # for file in file_chain:
59 
60  number_entries = 0
61  try:
62  number_entries = file_chain.GetEntries()
63  except AttributeError:
64  print('Could not load input file(s) %s.' % options.input_file)
65 
66  if number_entries == 0:
67  print('Data tree is empty or does not exist in file(s) %s. Exit.' % options.input_file)
68  sys.exit(0)
69 
70  # open the output root file
71  output_root_file = TFile(options.output_file, 'recreate')
72 
73  # create and draw histograms
74  gStyle.SetOptStat(0)
75  draw_exthits(file_chain)
76  draw_likelihoods(file_chain)
77 
78  # close the output file
79  output_root_file.Write()
80  output_root_file.Close()
81 
82 
83 def draw_exthits(file_chain):
84  """
85  Draw the ExtHit-related distributions.
86  """
87 
88  # NOTE: *.Draw() must precede *.GetListOfFunctions().Add() or the latter will be discarded!
89  detectorID = TH1F('DetectorID', 'Detector ID for ExtHits', 8, -0.5, 7.5)
90  file_chain.Draw('ExtHits.m_DetectorID&0x0F>>DetectorID', '')
91  detectorID.GetXaxis().SetTitle('0=undefined, 1=PXD, 2=SVD, 3=CDC, 4=TOP, 5=ARICH, 6=ECL, 7=KLM')
92  detectorID.GetListOfFunctions().Add(TNamed('Description', "0=undefined, 1=PXD, 2=SVD, 3=CDC, 4=TOP, 5=ARICH, 6=ECL, 7=KLM"))
93  detectorID.GetListOfFunctions().Add(TNamed('Check', "ECL > (KLM ~ TOP) >> ARICH; others ~ 0"))
94  detectorID.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
95  detectorID.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
96  detectorID.Write()
97 
98  tof = TH1F('TOF', 'Time of Flight for non-KLM ExtHits', 100, 0.0, 25.0)
99  file_chain.Draw('ExtHits.m_TOF>>TOF', '(ExtHits.m_DetectorID&0x0F)!=7')
100  tof.GetXaxis().SetTitle('t (ns)')
101  tof.GetListOfFunctions().Add(TNamed('Description', "Time of flight along track from IP to each ExtHit"))
102  tof.GetListOfFunctions().Add(TNamed('Check', "Dominant peak at ~4.5 ns, secondary peak at ~6 ns"))
103  tof.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
104  tof.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
105  tof.Write()
106 
107  tofKLM = TH1F('TOFKLM', 'Time of Flight for KLM ExtHits', 100, 0.0, 25.0)
108  file_chain.Draw('ExtHits.m_TOF>>TOFKLM', '(ExtHits.m_DetectorID&0x0F)==7')
109  tofKLM.GetXaxis().SetTitle('t (ns)')
110  tofKLM.GetListOfFunctions().Add(TNamed('Description', "Time of flight along track from IP to each ExtHit"))
111  tofKLM.GetListOfFunctions().Add(TNamed('Check', "Broad rising distribution at 7-11 ns, then falling to ~17 ns"))
112  tofKLM.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
113  tofKLM.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
114  tofKLM.Write()
115 
116  r = TH1F('r', 'r for non-KLM ExtHits', 40, 0.0, 200.0)
117  file_chain.Draw('ExtHits.getPosition().Perp()>>r', '(ExtHits.m_DetectorID&0x0F)!=7')
118  r.GetXaxis().SetTitle('r (cm)')
119  r.GetListOfFunctions().Add(TNamed('Description', "Radial position in transverse plane of each ExtHit"))
120  r.GetListOfFunctions().Add(TNamed('Check', "Sharp peak at 120 cm; broad peak at 140 cm"))
121  r.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
122  r.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
123  r.Write()
124 
125  z = TH1F('z', 'z for non-KLM ExtHits', 100, -200.0, 300.0)
126  file_chain.Draw('ExtHits.getPosition().Z()>>z', '(ExtHits.m_DetectorID&0x0F)!=7')
127  z.GetXaxis().SetTitle('z (cm)')
128  z.GetListOfFunctions().Add(TNamed('Description', "Axial position of each ExtHit"))
129  z.GetListOfFunctions().Add(TNamed('Check', "Broad peak centered at 0 cm; sharp dip at 0; secondary peak at 170"))
130  z.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
131  z.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
132  z.Write()
133 
134  rKLM = TH1F('rKLM', 'r for KLM ExtHits', 50, 100.0, 350.0)
135  file_chain.Draw('ExtHits.getPosition().Perp()>>rKLM', '(ExtHits.m_DetectorID&0x0F)==7')
136  rKLM.GetXaxis().SetTitle('r (cm)')
137  rKLM.GetListOfFunctions().Add(TNamed('Description', "Radial position in transverse plane of each ExtHit"))
138  rKLM.GetListOfFunctions().Add(TNamed('Check', "Low shoulder below 200 cm (EKLM); comb-like pattern above 200 cm (BKLM)"))
139  rKLM.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
140  rKLM.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
141  rKLM.Write()
142 
143  zKLM = TH1F('zKLM', 'z for KLM ExtHits', 140, -300.0, 400.0)
144  file_chain.Draw('ExtHits.getPosition().Z()>>zKLM', '(ExtHits.m_DetectorID&0x0F)==7')
145  zKLM.GetXaxis().SetTitle('z (cm)')
146  zKLM.GetListOfFunctions().Add(TNamed('Description', "Axial position of each ExtHit"))
147  zKLM.GetListOfFunctions().Add(TNamed('Check', "Comb-like pattern at z<-180 and z>280 cm (EKLM); broad peak at 0 (BKLM)"))
148  zKLM.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
149  zKLM.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
150  zKLM.Write()
151 
152  xy = TH2F('xy', 'y vs x for non-KLM ExtHits', 100, -200.0, 200.0, 100, -200.0, 200.0)
153  file_chain.Draw('ExtHits.getPosition().Y():ExtHits.getPosition().X()>>xy', '(ExtHits.m_DetectorID&0x0F)!=7')
154  xy.GetXaxis().SetTitle('x (cm)')
155  xy.GetYaxis().SetTitle('y (cm)')
156  xy.GetListOfFunctions().Add(TNamed('Description', "Position projected into transverse plane of each ExtHit"))
157  xy.GetListOfFunctions().Add(TNamed('Check', "Cylindrical, with most hits in ECL and TOP"))
158  xy.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
159  xy.GetListOfFunctions().Add(TNamed('MetaOptions', "box, expert"))
160  xy.Write()
161 
162  xz = TH2F('xz', 'x vs z for non-KLM ExtHits', 125, -200.0, 300.0, 100, -200.0, 200.0)
163  file_chain.Draw('ExtHits.getPosition().X():ExtHits.getPosition().Z()>>xz', '(ExtHits.m_DetectorID&0x0F)!=7')
164  xz.GetXaxis().SetTitle('z (cm)')
165  xz.GetYaxis().SetTitle('x (cm)')
166  xz.GetListOfFunctions().Add(TNamed('Description', "Position projected into x-z plane of each ExtHit"))
167  xz.GetListOfFunctions().Add(TNamed('Check', " "))
168  xz.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
169  xz.GetListOfFunctions().Add(TNamed('MetaOptions', "box, expert"))
170  xz.Write()
171 
172  yz = TH2F('yz', 'y vs z for non-KLM ExtHits', 125, -200.0, 300.0, 100, -200.0, 200.0)
173  file_chain.Draw('ExtHits.getPosition().Y():ExtHits.getPosition().Z()>>yz', '(ExtHits.m_DetectorID&0x0F)!=7')
174  yz.GetXaxis().SetTitle('z (cm)')
175  yz.GetYaxis().SetTitle('y (cm)')
176  yz.GetListOfFunctions().Add(TNamed('Description', "Position projected into y-z plane of each ExtHit"))
177  yz.GetListOfFunctions().Add(TNamed('Check', " "))
178  yz.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
179  yz.GetListOfFunctions().Add(TNamed('MetaOptions', "box, expert"))
180  yz.Write()
181 
182  xyKLM = TH2F('xyKLM', 'y vs x for KLM ExtHits', 140, -350.0, 350.0, 140, -350.0, 350.0)
183  file_chain.Draw('ExtHits.getPosition().Y():ExtHits.getPosition().X()>>xyKLM', '(ExtHits.m_DetectorID&0x0F)==7')
184  xyKLM.GetXaxis().SetTitle('x (cm)')
185  xyKLM.GetYaxis().SetTitle('y (cm)')
186  xyKLM.GetListOfFunctions().Add(TNamed('Description', "Position projected into transverse plane of each ExtHit"))
187  xyKLM.GetListOfFunctions().Add(TNamed('Check', "Octagonal (BKLM) and cylindrical (EKLM)"))
188  xyKLM.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
189  xyKLM.GetListOfFunctions().Add(TNamed('MetaOptions', "box, expert"))
190  xyKLM.Write()
191 
192  xzKLM = TH2F('xzKLM', 'x vs z for KLM ExtHits', 140, -300.0, 400.0, 140, -350.0, 350.0)
193  file_chain.Draw('ExtHits.getPosition().X():ExtHits.getPosition().Z()>>xzKLM', '(ExtHits.m_DetectorID&0x0F)==7')
194  xzKLM.GetXaxis().SetTitle('z (cm)')
195  xzKLM.GetYaxis().SetTitle('x (cm)')
196  xzKLM.GetListOfFunctions().Add(TNamed('Description', "Position projected into x-z plane of each ExtHit"))
197  xzKLM.GetListOfFunctions().Add(TNamed('Check', " "))
198  xzKLM.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
199  xzKLM.GetListOfFunctions().Add(TNamed('MetaOptions', "box, expert"))
200  xzKLM.Write()
201 
202  yzKLM = TH2F('yzKLM', 'y vs z for KLM ExtHits', 140, -300.0, 400.0, 140, -350.0, 350.0)
203  file_chain.Draw('ExtHits.getPosition().Y():ExtHits.getPosition().Z()>>yzKLM', '(ExtHits.m_DetectorID&0x0F)==7')
204  yzKLM.GetXaxis().SetTitle('z (cm)')
205  yzKLM.GetYaxis().SetTitle('y (cm)')
206  yzKLM.GetListOfFunctions().Add(TNamed('Description', "Position projected into y-z plane of each ExtHit"))
207  yzKLM.GetListOfFunctions().Add(TNamed('Check', " "))
208  yzKLM.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
209  yzKLM.GetListOfFunctions().Add(TNamed('MetaOptions', "box, expert"))
210  yzKLM.Write()
211 
212 
213 def draw_likelihoods(file_chain):
214  """
215  Draw the Muid likelihood-based distributions.
216  """
217 
218  outcome = TH1F('Outcome', 'Outcome', 5, -0.5, 4.5)
219 
220  llBins = 120
221  llMin = -60.0
222  llMax = 0.0
223  llMu_mu = TH1F('LLmu-muons', 'logL(mu) for true muons', llBins, llMin, llMax)
224  llMu_pi = TH1F('LLmu-pions', 'logL(mu) for true pions', llBins, llMin, llMax)
225  llPi_mu = TH1F('LLpi-muons', 'logL(pi) for true muons', llBins, llMin, llMax)
226  llPi_pi = TH1F('LLpi-pions', 'logL(pi) for true pions', llBins, llMin, llMax)
227  llDiff_mu = TH1F('LLdiff-muons', 'logL(mu)-logL(pi) for true muons', llBins * 2, llMin, -llMin)
228  llDiff_pi = TH1F('LLdiff-pions', 'logL(mu)-logL(pi) for true pions', llBins * 2, llMin, -llMin)
229 
230  layerBins = 15
231  layerMin = -0.5
232  layerMax = 14.5
233  barrelLayer_mu = TH1F('BarrelLayer-muons', 'Last barrel layer with hit for true muons', layerBins, layerMin, layerMax)
234  barrelLayer_pi = TH1F('BarrelLayer-pions', 'Last barrel layer with hit for true pions', layerBins, layerMin, layerMax)
235  endcapLayer_mu = TH1F('EndcapLayer-muons', 'Last endcap layer with hit for true muons', layerBins, layerMin, layerMax)
236  endcapLayer_pi = TH1F('EndcapLayer-pions', 'Last endcap layer with hit for true pions', layerBins, layerMin, layerMax)
237 
238  diffLayerBins = 20
239  diffLayerMin = -0.5
240  diffLayerMax = 19.5
241  diffLayer_mu = TH1F('DiffLayer-muons', 'Range difference for true muons', diffLayerBins, diffLayerMin, diffLayerMax)
242  diffLayer_pi = TH1F('DiffLayer-pions', 'Range difference for true pions', diffLayerBins, diffLayerMin, diffLayerMax)
243 
244  ndofBins = 32
245  ndofMin = -0.5
246  ndofMax = 31.5
247  ndof_mu = TH1F('Ndof-muons', 'Degrees of freedom for true muons', ndofBins, ndofMin, ndofMax)
248  ndof_pi = TH1F('Ndof-pions', 'Degrees of freedom for true pions', ndofBins, ndofMin, ndofMax)
249 
250  rchisqBins = 50
251  rchisqMin = 0.0
252  rchisqMax = 10.0
253  rchisq_mu = TH1F('rchisq-muons', 'Reduced chi-squared for true muons', rchisqBins, rchisqMin, rchisqMax)
254  rchisq_pi = TH1F('rchisq-pions', 'Reduced chi-squared for true pions', rchisqBins, rchisqMin, rchisqMax)
255 
256  momentumBins = 25 # for efficiency
257  momentumBins2 = 5 # for fake rate (low statistics)
258  momentumMin = 0.0
259  momentumMax = 5.0
260  efficiency_momentum = TH1F('Eff-momentum', 'Muon efficiency vs momentum', momentumBins, momentumMin, momentumMax)
261  efficiency_momentum_denom = TH1F('Eff-momentum-denom', 'Muon efficiency vs momentum', momentumBins, momentumMin, momentumMax)
262  fakerate_momentum = TH1F('FakeRate-momentum', 'Pion fake rate vs momentum', momentumBins2, momentumMin, momentumMax)
263  fakerate_momentum_denom = TH1F('FakeRate-momentum-denom', 'Pion fake rate vs momentum', momentumBins2, momentumMin, momentumMax)
264 
265  thetaBins = 35 # for efficiency
266  thetaBins2 = 7 # for fake rate (low statistics)
267  thetaMin = 10.0
268  thetaMax = 150.0
269  efficiency_theta = TH1F('Eff-theta', 'Muon efficiency vs theta', thetaBins, thetaMin, thetaMax)
270  efficiency_theta_denom = TH1F('Eff-theta-denom', 'Muon efficiency vs theta', thetaBins, thetaMin, thetaMax)
271  fakerate_theta = TH1F('FakeRate-theta', 'Pion fake rate vs theta', thetaBins2, thetaMin, thetaMax)
272  fakerate_theta_denom = TH1F('FakeRate-theta-denom', 'Pion fake rate vs theta', thetaBins2, thetaMin, thetaMax)
273 
274  phiBins = 36 # for efficiency
275  phiBins2 = 8 # for fake rate (low statistics)
276  phiMin = 0.0
277  phiMax = 360.0
278  efficiency_phi = TH1F('Eff-phi', 'Muon efficiency vs phi', phiBins, phiMin, phiMax)
279  efficiency_phi_denom = TH1F('Eff-phi-denom', 'Muon efficiency vs phi', phiBins, phiMin, phiMax)
280  fakerate_phi = TH1F('FakeRate-phi', 'Pion fake rate vs phi', phiBins2, phiMin, phiMax)
281  fakerate_phi_denom = TH1F('FakeRate-phi-denom', 'Pion fake rate vs phi', phiBins2, phiMin, phiMax)
282 
283  for entry in file_chain:
284  mcps = entry.MCParticles
285  momentum = mcps[0].getMomentum()
286  muids = entry.KLMMuidLikelihoods
287  for i in range(muids.GetEntriesFast()):
288  outcome.Fill(muids[i].getOutcome())
289  if muids[i].getOutcome() > 0:
290  llMu = muids[i].getLogL_mu()
291  llPi = muids[i].getLogL_pi()
292  blayer = muids[i].getBarrelHitLayer()
293  elayer = muids[i].getEndcapHitLayer()
294  diffLayer = muids[i].getExtLayer() - muids[i].getHitLayer()
295  ndof = muids[i].getDegreesOfFreedom()
296  chisq = muids[i].getChiSquared()
297  rchisq = -1.0
298  if ndof > 0:
299  rchisq = chisq / ndof
300  p = momentum.Mag()
301  theta = momentum.Theta() * 180.0 / np.pi
302  phi = momentum.Phi() * 180.0 / np.pi
303  if phi < 0.0:
304  phi = phi + 360.0
305  if 'muon' in file_chain.GetCurrentFile().GetName():
306  llMu_mu.Fill(llMu)
307  llPi_mu.Fill(llPi)
308  llDiff_mu.Fill(llMu - llPi)
309  barrelLayer_mu.Fill(blayer)
310  endcapLayer_mu.Fill(elayer)
311  diffLayer_mu.Fill(diffLayer)
312  ndof_mu.Fill(ndof)
313  rchisq_mu.Fill(rchisq)
314  efficiency_momentum_denom.Fill(p)
315  efficiency_theta_denom.Fill(theta)
316  efficiency_phi_denom.Fill(phi)
317  if llMu > llPi:
318  efficiency_momentum.Fill(p)
319  efficiency_theta.Fill(theta)
320  efficiency_phi.Fill(phi)
321  else:
322  llMu_pi.Fill(llMu)
323  llPi_pi.Fill(llPi)
324  llDiff_pi.Fill(llMu - llPi)
325  barrelLayer_pi.Fill(blayer)
326  endcapLayer_pi.Fill(elayer)
327  diffLayer_pi.Fill(diffLayer)
328  ndof_pi.Fill(ndof)
329  rchisq_pi.Fill(rchisq)
330  fakerate_momentum_denom.Fill(p)
331  fakerate_theta_denom.Fill(theta)
332  fakerate_phi_denom.Fill(phi)
333  if llMu > llPi:
334  fakerate_momentum.Fill(p)
335  fakerate_theta.Fill(theta)
336  fakerate_phi.Fill(phi)
337 
338  for j in range(efficiency_momentum_denom.GetNbinsX()):
339  num = efficiency_momentum.GetBinContent(j + 1)
340  denom = efficiency_momentum_denom.GetBinContent(j + 1)
341  if denom > 0:
342  efficiency_momentum.SetBinContent(j + 1, num / denom)
343  efficiency_momentum.SetBinError(j + 1, math.sqrt(num * (denom - num) / (denom * denom * denom)))
344  for j in range(efficiency_theta_denom.GetNbinsX()):
345  num = efficiency_theta.GetBinContent(j + 1)
346  denom = efficiency_theta_denom.GetBinContent(j + 1)
347  if denom > 0:
348  efficiency_theta.SetBinContent(j + 1, num / denom)
349  efficiency_theta.SetBinError(j + 1, math.sqrt(num * (denom - num) / (denom * denom * denom)))
350  for j in range(efficiency_phi_denom.GetNbinsX()):
351  num = efficiency_phi.GetBinContent(j + 1)
352  denom = efficiency_phi_denom.GetBinContent(j + 1)
353  if denom > 0:
354  efficiency_phi.SetBinContent(j + 1, num / denom)
355  efficiency_phi.SetBinError(j + 1, math.sqrt(num * (denom - num) / (denom * denom * denom)))
356  for j in range(fakerate_momentum_denom.GetNbinsX()):
357  num = fakerate_momentum.GetBinContent(j + 1)
358  denom = fakerate_momentum_denom.GetBinContent(j + 1)
359  if denom > 0:
360  fakerate_momentum.SetBinContent(j + 1, num / denom)
361  fakerate_momentum.SetBinError(j + 1, math.sqrt(num * (denom - num) / (denom * denom * denom)))
362  for j in range(fakerate_theta_denom.GetNbinsX()):
363  num = fakerate_theta.GetBinContent(j + 1)
364  denom = fakerate_theta_denom.GetBinContent(j + 1)
365  if denom > 0:
366  fakerate_theta.SetBinContent(j + 1, num / denom)
367  fakerate_theta.SetBinError(j + 1, math.sqrt(num * (denom - num) / (denom * denom * denom)))
368  for j in range(fakerate_phi_denom.GetNbinsX()):
369  num = fakerate_phi.GetBinContent(j + 1)
370  denom = fakerate_phi_denom.GetBinContent(j + 1)
371  if denom > 0:
372  fakerate_phi.SetBinContent(j + 1, num / denom)
373  fakerate_phi.SetBinError(j + 1, math.sqrt(num * (denom - num) / (denom * denom * denom)))
374 
375  # NOTE: *.Fill() must precede *.GetListOfFunctions().Add() or the latter will be discarded!
376  outcome.GetListOfFunctions().Add(TNamed('Description', "0=not in KLM, 1/2=barrel/endcap stop, 3/4=barrel/endcap exit"))
377  outcome.GetListOfFunctions().Add(TNamed('Check', "Peak at 3; rest ~ flat"))
378  outcome.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
379  outcome.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
380  outcome.SetMinimum(0.0)
381 
382  llMu_mu.GetListOfFunctions().Add(TNamed('Description', "Log-likelihood of muon hypothesis for true muons"))
383  llMu_mu.GetListOfFunctions().Add(TNamed('Check', "Sharp peak at 0 with tail on negative side"))
384  llMu_mu.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
385  llMu_mu.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
386  llMu_pi.GetListOfFunctions().Add(TNamed('Description', "Log-likelihood of muon hypothesis for true pions"))
387  llMu_pi.GetListOfFunctions().Add(TNamed('Check', "Sharp double peak at ~52 with tail up to 0"))
388  llMu_pi.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
389  llMu_pi.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
390  llPi_mu.GetListOfFunctions().Add(TNamed('Description', "Log-likelihood of pion hypothesis for true muons"))
391  llPi_mu.GetListOfFunctions().Add(TNamed('Check', "Sharp peak at ~-42 with tail up to 0"))
392  llPi_mu.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
393  llPi_mu.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
394  llPi_pi.GetListOfFunctions().Add(TNamed('Description', "Log-likelihood of pion hypothesis for true pions"))
395  llPi_pi.GetListOfFunctions().Add(TNamed('Check', "Sharp peak at 0 with tail on negative side"))
396  llPi_pi.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
397  llPi_pi.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
398  llDiff_mu.GetListOfFunctions().Add(TNamed('Description', "Log-likelihood difference of muon-pion hypotheses for true muons"))
399  llDiff_mu.GetListOfFunctions().Add(TNamed('Check', "Sharp peak at ~+42 with tail on lower side"))
400  llDiff_mu.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
401  llDiff_mu.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
402  llDiff_pi.GetListOfFunctions().Add(TNamed('Description', "Log-likelihood difference of muon-pion hypotheses for true pions"))
403  llDiff_pi.GetListOfFunctions().Add(TNamed('Check', "Sharp peak at ~-52 with tail on higher side"))
404  llDiff_pi.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
405  llDiff_pi.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
406 
407  barrelLayer_mu.GetListOfFunctions().Add(TNamed('Description', "Outermost barrel layer with matching hit for true muons"))
408  barrelLayer_mu.GetListOfFunctions().Add(TNamed('Check', "Sharp peak at layer 14"))
409  barrelLayer_mu.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
410  barrelLayer_mu.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
411  barrelLayer_pi.GetListOfFunctions().Add(TNamed('Description', "Outermost barrel layer with matching hit for true pions"))
412  barrelLayer_pi.GetListOfFunctions().Add(TNamed('Check', "Peak at layers 0-1 with tail above that"))
413  barrelLayer_pi.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
414  barrelLayer_pi.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
415  endcapLayer_mu.GetListOfFunctions().Add(TNamed('Description', "Outermost endcap layer with matching hit for true muons"))
416  endcapLayer_mu.GetListOfFunctions().Add(TNamed('Check', "Sharp peaks at layers 11 (backward) and 13 (forward)"))
417  endcapLayer_mu.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
418  endcapLayer_mu.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
419  endcapLayer_pi.GetListOfFunctions().Add(TNamed('Description', "Outermost endcap layer with matching hit for true pions"))
420  endcapLayer_pi.GetListOfFunctions().Add(TNamed('Check', "Peak at layers 0-1 with tail above that"))
421  endcapLayer_pi.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
422  endcapLayer_pi.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
423 
424  diffLayer_mu.GetListOfFunctions().Add(TNamed('Description', "Outermost extrapolated-vs-hit layer difference for true muons"))
425  diffLayer_mu.GetListOfFunctions().Add(TNamed('Check', "Sharp peak at 0"))
426  diffLayer_mu.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
427  diffLayer_mu.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
428  diffLayer_mu.GetXaxis().SetTitle('Last Ext Layer - Last Hit Layer')
429  diffLayer_pi.GetListOfFunctions().Add(TNamed('Description', "Outermost extrapolated-vs-hit layer difference for true pions"))
430  diffLayer_pi.GetListOfFunctions().Add(TNamed('Check', "Sharp peak at 15 with tail below that"))
431  diffLayer_pi.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
432  diffLayer_pi.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
433  diffLayer_pi.GetXaxis().SetTitle('Last Ext Layer - Last Hit Layer')
434 
435  ndof_mu.GetListOfFunctions().Add(TNamed('Description', "(Number of matching-hit layers) x 2 for true muons"))
436  ndof_mu.GetListOfFunctions().Add(TNamed('Check', "Even-only entries; peak at 30 with tail below"))
437  ndof_mu.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
438  ndof_mu.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
439  ndof_mu.GetXaxis().SetTitle('#dof')
440  ndof_pi.GetListOfFunctions().Add(TNamed('Description', "(Number of matching-hit layers) x 2 for true pions"))
441  ndof_pi.GetListOfFunctions().Add(TNamed('Check', "Even-only entries; sharp peak at 0 with tail above"))
442  ndof_pi.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
443  ndof_pi.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
444  ndof_pi.GetXaxis().SetTitle('#dof')
445 
446  rchisq_mu.GetListOfFunctions().Add(TNamed('Description', "Transverse-deviation chi-squared/#dof for true muons"))
447  rchisq_mu.GetListOfFunctions().Add(TNamed('Check', "Narrow peak at 1"))
448  rchisq_mu.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
449  rchisq_mu.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
450  rchisq_mu.GetXaxis().SetTitle('#chi^{2} / #ndof')
451  rchisq_pi.GetListOfFunctions().Add(TNamed('Description', "Transverse-deviation chi-squared/#dof for true pions"))
452  rchisq_pi.GetListOfFunctions().Add(TNamed('Check', "Broad peak at ~1"))
453  rchisq_pi.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
454  rchisq_pi.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=0.50,pvalue-error=0.10'))
455  rchisq_pi.GetXaxis().SetTitle('#chi^{2} / #ndof')
456 
457  efficiency_momentum.GetListOfFunctions().Add(TNamed('Description', "Muon ID efficiency vs p for logL(mu) > logL(pi)"))
458  efficiency_momentum.GetListOfFunctions().Add(TNamed('Check', "90-95% for p > 0.6 GeV/c"))
459  efficiency_momentum.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
460  efficiency_momentum.GetListOfFunctions().Add(TNamed('MetaOptions', 'pvalue-warn=0.50,pvalue-error=0.10'))
461  efficiency_momentum.GetXaxis().SetTitle('p (GeV/#it{c})')
462  efficiency_momentum.SetMinimum(0.0)
463  fakerate_momentum.GetListOfFunctions().Add(TNamed('Description', "Pion fake rate vs p for logL(mu) > logL(pi)"))
464  fakerate_momentum.GetListOfFunctions().Add(TNamed('Check', "Peak of 11% below 1 GeV/c, 2-6% above 1 GeV/c"))
465  fakerate_momentum.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
466  fakerate_momentum.GetListOfFunctions().Add(TNamed('MetaOptions', 'pvalue-warn=0.50,pvalue-error=0.10'))
467  fakerate_momentum.GetXaxis().SetTitle('p (GeV/#it{c})')
468  fakerate_momentum.SetMinimum(0.0)
469 
470  efficiency_theta.GetListOfFunctions().Add(TNamed('Description', "Muon ID efficiency vs theta for logL(mu) > logL(pi)"))
471  efficiency_theta.GetListOfFunctions().Add(TNamed('Check', "Somewhat flat (80-100%) in KLM range (17 < theta < 145 degrees)"))
472  efficiency_theta.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
473  efficiency_theta.GetListOfFunctions().Add(TNamed('MetaOptions', 'pvalue-warn=0.50,pvalue-error=0.10'))
474  efficiency_theta.GetXaxis().SetTitle('#theta (deg)')
475  efficiency_theta.SetMinimum(0.0)
476  fakerate_theta.GetListOfFunctions().Add(TNamed('Description', "Pion fake rate vs theta for logL(mu) > logL(pi)"))
477  fakerate_theta.GetListOfFunctions().Add(TNamed('Check', "Roughly flat (3-7%) in KLM range (17 < theta < 145 degrees)"))
478  fakerate_theta.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
479  fakerate_theta.GetListOfFunctions().Add(TNamed('MetaOptions', 'pvalue-warn=0.50,pvalue-error=0.10'))
480  fakerate_theta.GetXaxis().SetTitle('#theta (deg)')
481  fakerate_theta.SetMinimum(0.0)
482 
483  efficiency_phi.GetListOfFunctions().Add(TNamed('Description', "Muon identification efficiency vs phi for logL(mu) > logL(pi)"))
484  efficiency_phi.GetListOfFunctions().Add(TNamed('Check', "Somewhat flat distribution (85-100%)"))
485  efficiency_phi.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
486  efficiency_phi.GetListOfFunctions().Add(TNamed('MetaOptions', 'pvalue-warn=0.50,pvalue-error=0.10'))
487  efficiency_phi.GetXaxis().SetTitle('#phi (deg)')
488  efficiency_phi.SetMinimum(0.0)
489  fakerate_phi.GetListOfFunctions().Add(TNamed('Description', "Pion fake rate vs phi for logL(mu) > logL(pi)"))
490  fakerate_phi.GetListOfFunctions().Add(TNamed('Check', "Roughly flat distribution (2-7%)"))
491  fakerate_phi.GetListOfFunctions().Add(TNamed('Contact', "piilonen@vt.edu"))
492  fakerate_phi.GetListOfFunctions().Add(TNamed('MetaOptions', 'pvalue-warn=0.50,pvalue-error=0.10'))
493  fakerate_phi.GetXaxis().SetTitle('#phi (deg)')
494  fakerate_phi.SetMinimum(0.0)
495 
496  outcome.Write()
497  llMu_mu.Write()
498  llMu_pi.Write()
499  llPi_mu.Write()
500  llPi_pi.Write()
501  llDiff_mu.Write()
502  llDiff_pi.Write()
503  barrelLayer_mu.Write()
504  barrelLayer_pi.Write()
505  endcapLayer_mu.Write()
506  endcapLayer_pi.Write()
507  diffLayer_mu.Write()
508  diffLayer_pi.Write()
509  ndof_mu.Write()
510  ndof_pi.Write()
511  rchisq_mu.Write()
512  rchisq_pi.Write()
513  efficiency_momentum.Write()
514  efficiency_theta.Write()
515  efficiency_phi.Write()
516  fakerate_momentum.Write()
517  fakerate_theta.Write()
518  fakerate_phi.Write()
519 
520 
523 if __name__ == '__main__':
524  main()
main
int main(int argc, char **argv)
Run all tests.
Definition: test_main.cc:77