Belle II Software  release-05-01-25
BKLMMuonPlot.py
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 
4 
13 
14 """
15 <header>
16  <input>muon-KLMValidation.root</input>
17  <contact>martina.laurenza@roma3.infn.it</contact>
18  <description>Create validation plots for BKLM</description>
19 </header>
20 """
21 
22 import ROOT
23 ROOT.gROOT.SetBatch(ROOT.kTRUE) # noqa
24 
25 from ROOT import TFile, TChain, TH1F, TH2F, TNamed, gStyle, PyConfig
26 PyConfig.IgnoreCommandLineOptions = True # noqa
27 
28 import sys
29 import argparse
30 
31 # contact person information
32 # is added to the plot descriptions
33 CONTACT_PERSON = {'Name': 'Martina Laurenza',
34  'Email': 'martina.laurenza@roma3.infn.it'}
35 
36 
37 def main():
38  """
39  Create validation plots for BKLM.
40  """
41  print('Creating the BKLM validation plots...')
42 
43  parser = argparse.ArgumentParser()
44 
45  parser.add_argument('-i', '--input-file', dest='input_file',
46  default='../muon-KLMValidation.root',
47  help='Root file with Ext/Muid/BKLM/EKLM validation data.'
48  )
49  parser.add_argument('-o', '--output-file', dest='output_file',
50  default='BKLMMuon.root',
51  help='Root file with BKLM validation histograms.')
52 
53  args = parser.parse_args()
54 
55  # load chain of input files
56  file_chain = TChain('tree')
57  file_chain.Add(args.input_file)
58 
59  number_entries = 0
60  try:
61  number_entries = file_chain.GetEntries()
62  except AttributeError:
63  print('Could not load input file(s) %s.' % args.input_file)
64 
65  if number_entries == 0:
66  print('Data tree is empty or does not exist in file(s) %s. Exit.' % args.input_file)
67  sys.exit(0)
68 
69  # open the output root file
70  output_root_file = TFile(args.output_file, 'recreate')
71 
72  # create and draw histograms
73  gStyle.SetOptStat(0)
74  draw_bklmhists(file_chain)
75 
76  # close the output file
77  output_root_file.Write()
78  output_root_file.Close()
79 
80  print('BKLM validation plots created.')
81 
82 
83 # NOTE: *.Draw() must precede *.GetListOfFunctions().Add() or the latter will be discarded!
84 def draw_bklmhists(file_chain):
85  """
86  Draw the BKLMHit2d-related distributions.
87  """
88 
89  contact = 'Martina Laurenza (martina.laurenza@roma3.infn.it)'
90 
91  # Shifter plots
92 
93  inRPC = TH1F('InRPC', 'InRPC for BKLMHit2ds', 2, -0.5, 1.5)
94  file_chain.Draw('BKLMHit2ds.inRPC()>>InRPC', '')
95  inRPC.GetXaxis().SetTitle('0=scintillator 1=RPC')
96  inRPC.GetListOfFunctions().Add(TNamed('Description', 'Flag indicating if a muon hit is in scintillator (0) or RPC (1)'))
97  inRPC.GetListOfFunctions().Add(TNamed('Check', 'Mostly in RPC'))
98  inRPC.GetListOfFunctions().Add(TNamed('Contact', contact))
99  inRPC.GetListOfFunctions().Add(TNamed('MetaOptions', 'shifter,pvalue-warn=0.50,pvalue-error=0.10'))
100  inRPC.SetMinimum(0.0)
101  inRPC.Write()
102 
103  section = TH1F('Forward', 'Section for BKLMHit2ds', 2, -0.5, 1.5)
104  file_chain.Draw('BKLMHit2ds.getSection()>>Forward', '')
105  section.GetXaxis().SetTitle('0=backward 1=forward')
106  section.GetListOfFunctions().Add(TNamed('Description',
107  'Flag indicating if a muon hit is in backward (0) or forward (1) BKLM'))
108  section.GetListOfFunctions().Add(TNamed('Check', 'Somewhat more hits in the backward'))
109  section.GetListOfFunctions().Add(TNamed('Contact', contact))
110  section.GetListOfFunctions().Add(TNamed('MetaOptions', 'shifter,pvalue-warn=0.50,pvalue-error=0.10'))
111  section.SetMinimum(0.0)
112  section.Write()
113 
114  isOnTrack = TH1F('IsOnTrack', 'IsOnTrack for BKLMHit2ds', 2, -0.5, 1.5)
115  file_chain.Draw('BKLMHit2ds.isOnTrack()>>IsOnTrack', '')
116  isOnTrack.GetXaxis().SetTitle('0=not associated with Track 1=associated with Track')
117  isOnTrack.GetListOfFunctions().Add(TNamed('Description',
118  'Flag indicating if a muon hit is associated with a CDC Track by Muid'))
119  isOnTrack.GetListOfFunctions().Add(TNamed('Check', 'Mostly associated'))
120  isOnTrack.GetListOfFunctions().Add(TNamed('Contact', contact))
121  isOnTrack.GetListOfFunctions().Add(TNamed('MetaOptions', 'shifter,pvalue-warn=0.50,pvalue-error=0.10'))
122  isOnTrack.SetMinimum(0.0)
123  isOnTrack.Write()
124 
125  sector = TH1F('Sector', 'Sector for BKLMHit2ds', 10, -0.5, 9.5)
126  file_chain.Draw('BKLMHit2ds.getSector()>>Sector', '')
127  sector.GetXaxis().SetTitle('Sector #')
128  sector.GetListOfFunctions().Add(TNamed('Description', 'Sector number of muon hit'))
129  sector.GetListOfFunctions().Add(TNamed('Check', 'Roughly flat in sectors 1-8'))
130  sector.GetListOfFunctions().Add(TNamed('Contact', contact))
131  sector.GetListOfFunctions().Add(TNamed('MetaOptions', 'shifter,pvalue-warn=1.00,pvalue-error=0.01'))
132  sector.SetMinimum(0.0)
133  sector.Write()
134 
135  layer = TH1F('Layer', 'Layer for BKLMHit2ds', 16, -0.5, 15.5)
136  file_chain.Draw('BKLMHit2ds.getLayer()>>Layer', '')
137  layer.GetXaxis().SetTitle('Layer #')
138  layer.GetListOfFunctions().Add(TNamed('Description', 'Layer number of muon hit'))
139  layer.GetListOfFunctions().Add(TNamed('Check',
140  'First peak at layer 1 and second (higher) peak at layer 2, with tail above those'))
141  layer.GetListOfFunctions().Add(TNamed('Contact', contact))
142  layer.GetListOfFunctions().Add(TNamed('MetaOptions', 'shifter,pvalue-warn=1.0,pvalue-error=0.01'))
143  layer.SetMinimum(0.0)
144  layer.Write()
145 
146  phistrip = TH1F('PhiStrip', 'PhiStrip for BKLMHit2ds', 50, -0.5, 49.5)
147  file_chain.Draw('BKLMHit2ds.getPhiStripAve()>>PhiStrip', '')
148  phistrip.GetXaxis().SetTitle('Phi strip #')
149  phistrip.GetListOfFunctions().Add(TNamed('Description', 'Strip number in phi plane of muon hit'))
150  phistrip.GetListOfFunctions().Add(TNamed('Check',
151  'Roughly flat for 1-36 (all layers) and then for 37-48 (layers 6-15)'))
152  phistrip.GetListOfFunctions().Add(TNamed('Contact', contact))
153  phistrip.GetListOfFunctions().Add(TNamed('MetaOptions', 'shifter,pvalue-warn=1.00,pvalue-error=0.01'))
154  phistrip.SetMinimum(0.0)
155  phistrip.Write()
156 
157  zstrip = TH1F('ZStrip', 'ZStrip for BKLMHit2ds', 60, -0.5, 59.5)
158  file_chain.Draw('BKLMHit2ds.getZStripAve()>>ZStrip', '')
159  zstrip.GetXaxis().SetTitle('Z strip #')
160  zstrip.GetListOfFunctions().Add(TNamed('Description', 'Strip number in z plane of muon hit'))
161  zstrip.GetListOfFunctions().Add(TNamed('Check',
162  'Downward-sloping for 1-48 (all layers), shoulder for 49-54 (layers 1-2)'))
163  zstrip.GetListOfFunctions().Add(TNamed('Contact', contact))
164  zstrip.GetListOfFunctions().Add(TNamed('MetaOptions', 'shifter,pvalue-warn=1.00,pvalue-error=0.01'))
165  zstrip.SetMinimum(0.0)
166  zstrip.Write()
167 
168  RPC_bklmdigits_tres = ROOT.TH1F('RPC_bklmdigits_tres', 'KLM Digits time resolution in RPC', 250, -25, 25)
169  file_chain.Draw('KLMDigits.getTime()-KLMDigits.getMCTime()>>RPC_bklmdigits_tres',
170  'KLMDigits.getSubdetector()==1 && KLMDigits.getLayer() > 3')
171  RPC_bklmdigits_tres.SetXTitle('ns')
172  RPC_bklmdigits_tres.GetListOfFunctions().Add(TNamed('Description', 'KLMDigits Time resolution in RPC'))
173  RPC_bklmdigits_tres.GetListOfFunctions().Add(TNamed('Contact', contact))
174  RPC_bklmdigits_tres.GetListOfFunctions().Add(TNamed('MetaOptions', 'shifter,pvalue-warn=0.50,pvalue-error=0.10'))
175  RPC_bklmdigits_tres.Write()
176 
177  timeRPC = TH1F('TimeRPC', 'Hit time for BKLMHit2ds in RPCs', 100, -2.0, 2.0)
178  file_chain.Draw('BKLMHit2ds.getTime()>>TimeRPC', 'BKLMHit2ds.inRPC()==1')
179  timeRPC.GetXaxis().SetTitle('t (ns)')
180  timeRPC.GetListOfFunctions().Add(TNamed('Description', 'Time of muon hit in RPCs'))
181  timeRPC.GetListOfFunctions().Add(TNamed('Check', 'Narrow peak at 0 ns'))
182  timeRPC.GetListOfFunctions().Add(TNamed('Contact', contact))
183  timeRPC.GetListOfFunctions().Add(TNamed('MetaOptions', 'shifter,pvalue-warn=1.00,pvalue-error=0.01'))
184  timeRPC.Write()
185 
186  Sci_bklmdigits_tres = ROOT.TH1F('Sci_bklmdigits_tres', 'KLM Digits time resolution in Scintillators', 250, -25, 25)
187  file_chain.Draw('KLMDigits.getTime()-KLMDigits.getMCTime()>>Sci_bklmdigits_tres',
188  'KLMDigits.getSubdetector()==1 && KLMDigits.getLayer() <= 3')
189  Sci_bklmdigits_tres.SetXTitle('ns')
190  Sci_bklmdigits_tres.GetListOfFunctions().Add(TNamed('Description', 'KLMDigits Time resolution in Scintillators'))
191  Sci_bklmdigits_tres.GetListOfFunctions().Add(TNamed('Contact', contact))
192  Sci_bklmdigits_tres.GetListOfFunctions().Add(TNamed('MetaOptions', 'shifter,pvalue-warn=0.50,pvalue-error=0.10'))
193  Sci_bklmdigits_tres.Write()
194 
195  timeSci = TH1F('TimeSci', 'Hit time for BKLMHit2ds in scintillators', 100, -5.0, 15.0)
196  file_chain.Draw('BKLMHit2ds.getTime()>>TimeSci', 'BKLMHit2ds.inRPC()==0')
197  timeSci.GetXaxis().SetTitle('t (ns)')
198  timeSci.GetListOfFunctions().Add(TNamed('Description', 'Time of muon hit in scintillators'))
199  timeSci.GetListOfFunctions().Add(TNamed('Check',
200  'Broad peak mainly between 2 ns and 8 ns, with the mean around 3.5 ns'))
201  timeSci.GetListOfFunctions().Add(TNamed('Contact', contact))
202  timeSci.GetListOfFunctions().Add(TNamed('MetaOptions', 'shifter,pvalue-warn=1.00,pvalue-error=0.01'))
203  timeSci.Write()
204 
205  nPE = TH1F('nGenPE', 'Generated PE in BKLM', 50, 0.0, 100)
206  file_chain.Draw('KLMDigits.getNGeneratedPhotoelectrons()>>nGenPE', 'KLMDigits.getSubdetector()==1 && KLMDigits.m_Layer < 3')
207  nPE.GetXaxis().SetTitle('# generated PE')
208  nPE.GetListOfFunctions().Add(TNamed('Description', 'Number of generated photoelectrons in BKLM'))
209  nPE.GetListOfFunctions().Add(TNamed('Check', ''))
210  nPE.GetListOfFunctions().Add(TNamed('Contact', contact))
211  nPE.GetListOfFunctions().Add(TNamed('MetaOptions', ''))
212  nPE.Write()
213 
214  # Expert plots
215 
216  edep = TH1F('EnergyDeposit', 'Energy deposition for BKLMHit2ds', 50, 0.0, 25.0)
217  file_chain.Draw('BKLMHit2ds.getEnergyDeposit()*1000.0>>EnergyDeposit', '')
218  edep.GetXaxis().SetTitle('E (keV)')
219  edep.GetListOfFunctions().Add(TNamed('Description', 'dE/dx energy deposition of muon hit'))
220  edep.GetListOfFunctions().Add(TNamed('Check', 'Peak near 3 keV'))
221  edep.GetListOfFunctions().Add(TNamed('Contact', contact))
222  edep.GetListOfFunctions().Add(TNamed('MetaOptions', 'expert,pvalue-warn=1.00,pvalue-error=0.01'))
223  edep.Write()
224 
225 
226 # Entry point of this script: call the main() function
227 if __name__ == '__main__':
228  main()
main
int main(int argc, char **argv)
Run all tests.
Definition: test_main.cc:77