Belle II Software  release-05-01-25
laserPixelIlumination.py
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 
4 from basf2 import *
5 import sys
6 from ROOT import Belle2
7 from ROOT import TH2F, TFile
8 
9 # ------------------------------------------------------------------------
10 # example of making histograms of pixel hits coming from individual fibers
11 # needs as input the file produced by top/analysis/simLaserCalibration.py
12 # ------------------------------------------------------------------------
13 
14 inputFile = 'laserSimulation.root'
15 if not os.path.exists(inputFile):
16  B2ERROR(inputFile + ': file not found')
17  B2INFO('File can be generated with top/analysis/simLaserCalibration.py')
18  sys.exit()
19 
20 
21 class Histogrammer(Module):
22 
23  ''' A module to histogram pixel hits from individual fibers'''
24 
25 
26  hist = [TH2F('fiber' + str(k + 1), 'Pixel hit distribution from fiber No.' + str(k + 1),
27  64, 0.5, 64.5, 8, 0.5, 8.5) for k in range(9)]
28 
29  pixelID = 482
30 
31  propTime = TH2F('propTime', 'Photon propagation times for pixel ' + str(pixelID) +
32  ' vs. fiber number', 9, 0.5, 9.5, 200, 0.0, 1.0)
33  propTime.GetXaxis().SetTitle('fiber number')
34  propTime.GetYaxis().SetTitle('propagation time [ns]')
35 
36  def event(self):
37  ''' Event processor: fill histograms '''
38 
39  digits = Belle2.PyStoreArray('TOPDigits')
40  for digit in digits:
41  simhits = digit.getRelationsWith('TOPSimHits') # several simhits can contribute
42  for simhit in simhits:
43  photon = simhit.getRelated('TOPSimPhotons') # one or none possible
44  if photon:
45  k = int((photon.getEmissionPoint().X() + 22.5) / 5.5)
46  if k >= 0 and k < 9:
47  self.hist[k].Fill(digit.getPixelCol(), digit.getPixelRow())
48  else:
49  B2ERROR('wrong decoding of fiber number: ' + str(k + 1))
50  if digit.getPixelID() == self.pixelID:
51  t = photon.getDetectionTime() - photon.getEmissionTime()
52  self.propTime.Fill(k + 1, t)
53 
54  def terminate(self):
55  ''' Write histograms to file '''
56 
57  tfile = TFile('laserPixelIlumination.root', 'recreate')
58  for k in range(9):
59  self.hist[k].Write()
60  self.propTime.Write()
61  tfile.Close()
62 
63 
64 # Create path
65 main = create_path()
66 
67 # Input
68 roinput = register_module('RootInput')
69 roinput.param('inputFileName', inputFile)
70 main.add_module(roinput)
71 
72 # Histogrammer
73 main.add_module(Histogrammer())
74 
75 # Show progress of processing
76 progress = register_module('Progress')
77 main.add_module(progress)
78 
79 # Process events
80 process(main)
81 
82 # Print call statistics
83 print(statistics)
laserPixelIlumination.Histogrammer.hist
list hist
2D histograms
Definition: laserPixelIlumination.py:26
laserPixelIlumination.Histogrammer.terminate
def terminate(self)
Definition: laserPixelIlumination.py:54
laserPixelIlumination.Histogrammer.pixelID
int pixelID
pixel ID
Definition: laserPixelIlumination.py:29
laserPixelIlumination.Histogrammer.propTime
propTime
histogram
Definition: laserPixelIlumination.py:31
laserPixelIlumination.Histogrammer
Definition: laserPixelIlumination.py:21
Belle2::PyStoreArray
a (simplified) python wrapper for StoreArray.
Definition: PyStoreArray.h:58
laserPixelIlumination.Histogrammer.event
def event(self)
Definition: laserPixelIlumination.py:36