Belle II Software  release-05-01-25
EventDisplayer.py
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 
4 # Purpose:
5 # basf module to create BKLM event displays from BKLMHit2ds, ExtHits, and MuidHits.
6 #
7 
8 import basf2
9 import bklmDB
10 import rawdata
11 import math
12 import ctypes
13 import ROOT
14 from ROOT import Belle2, TH1F, TH2F, TCanvas, THistPainter, TPad
15 
16 # =========================================================================
17 #
18 # EventDisplayer
19 #
20 # =========================================================================
21 
22 
23 class EventDisplayer(basf2.Module):
24  """Draw BKLM event displays from BKLMHit2ds, ExtHits, and MuidHits."""
25 
26 
27  BKLM_ID = 0x07000000
28 
29  EKLM_ID = 0x08000000
30 
31  BKLM_STRIP_BIT = 0
32 
33  BKLM_PLANE_BIT = 6
34 
35  BKLM_LAYER_BIT = 7
36 
37  BKLM_SECTOR_BIT = 11
38 
39  BKLM_SECTION_BIT = 14
40 
41  BKLM_MAXSTRIP_BIT = 15
42 
43  BKLM_STRIP_MASK = 0x3f
44 
45  BKLM_PLANE_MASK = (1 << BKLM_PLANE_BIT)
46 
47  BKLM_LAYER_MASK = (15 << BKLM_LAYER_BIT)
48 
49  BKLM_SECTOR_MASK = (7 << BKLM_SECTOR_BIT)
50 
51  BKLM_SECTION_MASK = (1 << BKLM_SECTION_BIT)
52 
53  BKLM_MAXSTRIP_MASK = (63 << BKLM_MAXSTRIP_BIT)
54 
55  BKLM_MODULEID_MASK = (BKLM_SECTION_MASK | BKLM_SECTOR_MASK | BKLM_LAYER_MASK)
56 
57  def __init__(self, exp, run, eventPdfName, maxDisplays, minRPCHits, minMuidHits):
58  """Constructor
59 
60  Arguments:
61  exp (str): formatted experiment number
62  run (str): formatter run number
63  eventPdfName (str): path name of the output event-display PDF file
64  maxDisplays (int): max # of events displays to write
65  minRPCHits (int): min # of RPC BKLMHit2ds in any sector for event display
66  minMuidHits (int): min # of MuidHits in the event for event display
67  """
68  super().__init__()
69 
70  self.exp = exp
71 
72  self.run = run
73 
74  self.eventPdfName = eventPdfName
75 
76  self.maxDisplays = maxDisplays
77 
78  self.minRPCHits = minRPCHits
79 
80  self.minMuidHits = minMuidHits
81 
82  self.eventCounter = 0
83 
84  self.eventDisplays = 0
85 
86  self.lastTitle = ''
87 
88  def initialize(self):
89  """Handle job initialization: fill the mapping database, create histograms, open the event-display file"""
90 
91  expRun = 'e{0:02d}r{1}: '.format(int(self.exp), int(self.run))
92 
93  # Open the output PDF file for event displays
94 
95 
96  self.eventCanvas = ROOT.TCanvas("eventCanvas", self.eventPdfName, 3200, 1600)
97  title = '{0}['.format(self.eventPdfName)
98  self.eventCanvas.SaveAs(title)
99  self.eventCanvas.Clear()
100  self.eventCanvas.Divide(2, 1)
101 
102  # Create the boilerplate for the end- and side-views of the event display
103 
104 
105  self.cosine = [0, 0, 0, 0, 0, 0, 0, 0]
106 
107  self.sine = [0, 0, 0, 0, 0, 0, 0, 0]
108  for sector in range(0, 8):
109  phi = math.pi * sector / 4
110  self.cosine[sector] = math.cos(phi)
111  self.sine[sector] = math.sin(phi)
112 
113  self.hist_XY = ROOT.TH2F('XY', ' ;x;y', 10, -345.0, 345.0, 10, -345.0, 345.0)
114  self.hist_XY.SetStats(False)
115 
116  self.hist_ZY = ROOT.TH2F('ZY', ' ;z;y', 10, -345.0, 345.0, 10, -345.0, 345.0)
117  self.hist_ZY.SetStats(False)
118  # 300x300 cm^2 grid for each octant
119  u1 = 65
120  u2 = 365
121  u3 = 100
122  u4 = 400
123  u5 = 150
124 
125  self.hist_XYS = [0, 0, 0, 0, 0, 0, 0, 0]
126  self.hist_XYS[0] = ROOT.TH2F('XYS0', ' ;x;y', 10, +u3, +u4, 10, -u5, +u5)
127  self.hist_XYS[0].SetStats(False)
128  self.hist_XYS[1] = ROOT.TH2F('XYS1', ' ;x;y', 10, +u1, +u2, 10, +u1, +u2)
129  self.hist_XYS[1].SetStats(False)
130  self.hist_XYS[2] = ROOT.TH2F('XYS2', ' ;x;y', 10, -u5, +u5, 10, +u3, +u4)
131  self.hist_XYS[2].SetStats(False)
132  self.hist_XYS[3] = ROOT.TH2F('XYS3', ' ;x;y', 10, -u2, -u1, 10, +u1, +u2)
133  self.hist_XYS[3].SetStats(False)
134  self.hist_XYS[4] = ROOT.TH2F('XYS4', ' ;x;y', 10, -u4, -u3, 10, -u5, +u5)
135  self.hist_XYS[4].SetStats(False)
136  self.hist_XYS[5] = ROOT.TH2F('XYS5', ' ;x;y', 10, -u2, -u1, 10, -u2, -u1)
137  self.hist_XYS[5].SetStats(False)
138  self.hist_XYS[6] = ROOT.TH2F('XYS6', ' ;x;y', 10, -u5, +u5, 10, -u4, -u3)
139  self.hist_XYS[6].SetStats(False)
140  self.hist_XYS[7] = ROOT.TH2F('XYS7', ' ;x;y', 10, +u1, +u2, 10, -u2, -u1)
141  self.hist_XYS[7].SetStats(False)
142 
143  self.hist_ZYS = ROOT.TH2F('ZYS', ' ;z;y', 10, -150.0, 150.0, 10, 125.0, 425.0)
144  self.hist_ZYS.SetStats(False)
145  ROOT.gStyle.SetOptStat(10)
146 
147  self.bklmXY = []
148  r0 = 201.9 + 0.5 * 4.4 # cm
149  dr = 9.1 # cm
150  tan0 = math.tan(math.pi / 8.0)
151  g = ROOT.TGraph()
152  g.SetPoint(0, -200.0, 0.0)
153  g.SetPoint(1, +200.0, 0.0)
154  g.SetLineColor(19)
155  g.SetLineWidth(1)
156  g.SetLineStyle(3)
157  self.bklmXY.append(g)
158  g = ROOT.TGraph()
159  g.SetPoint(0, 0.0, -200.0)
160  g.SetPoint(1, 0.0, +200.0)
161  g.SetLineColor(19)
162  g.SetLineWidth(1)
163  g.SetLineStyle(3)
164  self.bklmXY.append(g)
165  g = ROOT.TGraph()
166  g.SetPoint(0, -5.0, 0.0)
167  g.SetPoint(1, +5.0, 0.0)
168  g.SetPoint(2, 0.0, 0.0)
169  g.SetPoint(3, 0.0, +5.0)
170  g.SetPoint(4, 0.0, -5.0)
171  g.SetLineColor(1)
172  g.SetLineWidth(1)
173  self.bklmXY.append(g)
174  for layer in range(0, 15):
175  r = r0 + layer * dr
176  x = r * tan0
177  g = ROOT.TGraph()
178  g.SetPoint(0, +r, -x)
179  g.SetPoint(1, +r, +x)
180  g.SetPoint(2, +x, +r)
181  g.SetPoint(3, -x, +r)
182  g.SetPoint(4, -r, +x)
183  g.SetPoint(5, -r, -x)
184  g.SetPoint(6, -x, -r)
185  g.SetPoint(7, +x, -r)
186  g.SetPoint(8, +r, -x)
187  if layer < 2:
188  g.SetLineColor(18)
189  else:
190  g.SetLineColor(17)
191  if (layer % 5) == 0:
192  g.SetLineStyle(3)
193  g.SetLineWidth(1)
194  self.bklmXY.append(g)
195 
196  self.bklmZY = []
197  rF = r0 + 14 * dr
198  x0 = r0 * tan0
199  z0 = 47.0 # cm
200  zL = 220.0 # cm
201  g = ROOT.TGraph()
202  g.SetPoint(0, -zL + z0 - 140.0, 0.0)
203  g.SetPoint(1, +zL + z0 + 70.0, 0.0)
204  g.SetLineColor(19)
205  g.SetLineWidth(1)
206  g.SetLineStyle(3)
207  self.bklmZY.append(g)
208  g = ROOT.TGraph()
209  g.SetPoint(0, 0.0, -315.0)
210  g.SetPoint(1, 0.0, +340.0)
211  g.SetLineColor(19)
212  g.SetLineWidth(1)
213  g.SetLineStyle(3)
214  self.bklmZY.append(g)
215  g = ROOT.TGraph()
216  g.SetPoint(0, -5.0, 0.0)
217  g.SetPoint(1, +5.0, 0.0)
218  g.SetPoint(2, 0.0, 0.0)
219  g.SetPoint(3, 0.0, +5.0)
220  g.SetPoint(4, 0.0, -5.0)
221  g.SetLineColor(1)
222  g.SetLineWidth(1)
223  self.bklmZY.append(g)
224  g = ROOT.TGraph()
225  g.SetPoint(0, -zL + z0, +x0)
226  g.SetPoint(1, -zL + z0, +r0)
227  g.SetLineColor(18)
228  g.SetLineWidth(1)
229  g.SetLineStyle(3)
230  self.bklmZY.append(g)
231  g = ROOT.TGraph()
232  g.SetPoint(0, -zL + z0, -x0)
233  g.SetPoint(1, -zL + z0, -r0)
234  g.SetLineColor(18)
235  g.SetLineWidth(1)
236  g.SetLineStyle(3)
237  self.bklmZY.append(g)
238  g = ROOT.TGraph()
239  g.SetPoint(0, +zL + z0, +x0)
240  g.SetPoint(1, +zL + z0, +r0)
241  g.SetLineColor(18)
242  g.SetLineWidth(1)
243  g.SetLineStyle(3)
244  self.bklmZY.append(g)
245  g = ROOT.TGraph()
246  g.SetPoint(0, +zL + z0, -x0)
247  g.SetPoint(1, +zL + z0, -r0)
248  g.SetLineColor(18)
249  g.SetLineWidth(1)
250  g.SetLineStyle(3)
251  self.bklmZY.append(g)
252  g = ROOT.TGraph()
253  g.SetPoint(0, -zL + z0, r0)
254  g.SetPoint(1, +zL + z0, r0)
255  g.SetPoint(2, +zL + z0, rF)
256  g.SetPoint(3, -zL + z0, rF)
257  g.SetPoint(4, -zL + z0, r0)
258  g.SetLineColor(18)
259  g.SetLineWidth(1)
260  self.bklmZY.append(g)
261  g = ROOT.TGraph()
262  g.SetPoint(0, -zL + z0, -r0)
263  g.SetPoint(1, +zL + z0, -r0)
264  g.SetPoint(2, +zL + z0, -rF)
265  g.SetPoint(3, -zL + z0, -rF)
266  g.SetPoint(4, -zL + z0, -r0)
267  g.SetLineColor(18)
268  g.SetLineWidth(1)
269  self.bklmZY.append(g)
270  g = ROOT.TGraph()
271  g.SetPoint(0, -zL + z0, -x0)
272  g.SetPoint(1, +zL + z0, -x0)
273  g.SetPoint(2, +zL + z0, +x0)
274  g.SetPoint(3, -zL + z0, +x0)
275  g.SetPoint(4, -zL + z0, -x0)
276  g.SetLineColor(18)
277  g.SetLineWidth(1)
278  self.bklmZY.append(g)
279 
280  self.bklmZYL = []
281  for layer in range(0, 15):
282  r = r0 + layer * dr
283  z0 = 47.0 # cm
284  zL = 220.0 # cm
285  g = ROOT.TGraph()
286  g.SetPoint(0, -zL + z0, r)
287  g.SetPoint(1, +zL + z0, r)
288  g.SetLineColor(19)
289  g.SetLineWidth(1)
290  self.bklmZYL.append(g)
291 
293 
294  self.sectorFBToDC = [11, 15, 2, 6, 10, 14, 3, 7, 9, 13, 0, 4, 8, 12, 1, 5]
295 
296  self.dcToSectorFB = [10, 14, 2, 6, 11, 15, 3, 7, 12, 8, 4, 0, 13, 9, 5, 1]
297 
298  self.t0Cal = 312
299 
300  self.t0Cal2d = 308
301 
302  self.ct0Cal = 455
303 
304  self.ct0Cal2d = 520
305 
306  self.t0RPC = [8, -14, -6, -14, -2, 10, 9, 13, 0, -10, -14, -20, 2, 6, 14, 11]
307 
308  self.ct0Scint = [-1, -33, -46, -33, -2, 32, 51, 32, 0, -32, -45, -33, -4, 34, 45, 27]
309 
310  def terminate(self):
311  """Handle job termination: close event-display file"""
312 
313  pdfNameLast = '{0}]'.format(self.eventPdfName)
314  self.eventCanvas.Print(pdfNameLast, self.lastTitle)
315  print('Goodbye')
316 
317  def beginRun(self):
318  """Handle begin of run: print diagnostic message"""
319  EventMetaData = Belle2.PyStoreObj('EventMetaData')
320  print('beginRun', EventMetaData.getRun())
321 
322  def endRun(self):
323  """Handle end of run: print diagnostic message"""
324  EventMetaData = Belle2.PyStoreObj('EventMetaData')
325  print('endRun', EventMetaData.getRun())
326 
327  def event(self):
328  """Process one event: (optionally) draw event display"""
329 
330  super().return_value(self.eventDisplays < self.maxDisplays)
331 
332  if self.eventDisplays >= self.maxDisplays:
333  return
334 
335  self.eventCounter += 1
336  EventMetaData = Belle2.PyStoreObj('EventMetaData')
337  event = EventMetaData.getEvent()
338  rawklms = Belle2.PyStoreArray('RawKLMs') # to determine if BKLMHit2d is prompt
339  hit2ds = Belle2.PyStoreArray('BKLMHit2ds')
340  exthits = Belle2.PyStoreArray('ExtHits')
341  muidhits = Belle2.PyStoreArray('MuidHits')
342 
343  rawFb = [[], [], [], [], [], [], [], [], [], [], [], [], [], [], [], []]
344  rawSector = [[], [], [], [], [], [], [], [], [], [], [], [], [], [], [], []]
345  rawLayer = [[], [], [], [], [], [], [], [], [], [], [], [], [], [], [], []]
346  rawPlane = [[], [], [], [], [], [], [], [], [], [], [], [], [], [], [], []]
347  rawStrip = [[], [], [], [], [], [], [], [], [], [], [], [], [], [], [], []]
348  rawCtime = [[], [], [], [], [], [], [], [], [], [], [], [], [], [], [], []]
349 
350  for copper in range(0, len(rawklms)):
351  rawklm = rawklms[copper]
352  nodeID = rawklm.GetNodeID(0) - self.BKLM_ID
353  if nodeID >= self.EKLM_ID - self.BKLM_ID:
354  nodeID = nodeID - (self.EKLM_ID - self.BKLM_ID) + 4
355  if (nodeID < 0) or (nodeID > 4): # skip EKLM nodes
356  continue
357  for finesse in range(0, 4):
358  dc = (finesse << 2) + copper
359  n = rawklm.GetDetectorNwords(0, finesse)
360  bufSlot = rawklm.GetDetectorBuffer(0, finesse)
361  if n > 0:
362  n = n >> 1
363  for j in range(0, n):
364  word0 = bufSlot[j * 2]
365  word1 = bufSlot[j * 2 + 1]
366  ctime = word0 & 0xffff
367  channel = (word0 >> 16) & 0x7f
368  axis = (word0 >> 23) & 0x01
369  lane = (word0 >> 24) & 0x1f # crate's slot number
370  flag = (word0 >> 30) & 0x03
371  fb = -1
372  sector = -1
373  layer = -1
374  plane = -1
375  strip = -1
376  electId = (channel << 12) | (axis << 11) | (lane << 6) | (finesse << 4) | nodeID
377  if electId in self.electIdToModuleId:
378  moduleId = self.electIdToModuleId[electId]
379  fb = (moduleId & self.BKLM_SECTION_MASK) >> self.BKLM_SECTION_BIT
380  sector = (moduleId & self.BKLM_SECTOR_MASK) >> self.BKLM_SECTOR_BIT
381  layer = (moduleId & self.BKLM_LAYER_MASK) >> self.BKLM_LAYER_BIT
382  plane = (moduleId & self.BKLM_PLANE_MASK) >> self.BKLM_PLANE_BIT
383  strip = (moduleId & self.BKLM_STRIP_MASK) >> self.BKLM_STRIP_BIT
384  rawFb[dc].append(fb)
385  rawSector[dc].append(sector)
386  rawLayer[dc].append(layer)
387  rawPlane[dc].append(plane)
388  rawStrip[dc].append(strip)
389  rawCtime[dc].append(ctime)
390 
391  tCal2d = []
392  for hit2d in hit2ds:
393  tCal2d.append(hit2d.getTime())
394 
395  # Process the ExtHits for event display
396 
397  xold = 0
398  yold = 0
399  zold = 0
400  sold = 0
401  epositions = []
402  sumx = 0
403  sumy = 0
404  sumz = 0
405  sumn = 0
406  for exthit in exthits:
407  sector = (exthit.getCopyID() & self.BKLM_SECTOR_MASK) >> self.BKLM_SECTOR_BIT
408  extPosition = exthit.getPosition()
409  x = extPosition[0]
410  y = extPosition[1]
411  z = extPosition[2]
412  dx = x - xold
413  dy = y - yold
414  dz = z - zold
415  if (dx * dx + dy * dy + dz * dz > 36) and (sumn > 0):
416  eposition = [sumx / sumn, sumy / sumn, sumz / sumn, sold]
417  epositions.append(eposition)
418  sumx = 0
419  sumy = 0
420  sumz = 0
421  sumn = 0
422  else:
423  if sumn == 0:
424  xold = x
425  yold = y
426  zold = z
427  sold = sector
428  sumx = sumx + x
429  sumy = sumy + y
430  sumz = sumz + z
431  sumn = sumn + 1
432  if sumn > 0:
433  eposition = [sumx / sumn, sumy / sumn, sumz / sumn, sold]
434  epositions.append(eposition)
435  extXYGraph = ROOT.TGraph()
436  extXYGraph.SetMarkerColor(30)
437  extXYGraph.SetMarkerSize(2.25)
438  extXYGraph.SetMarkerStyle(21)
439  extZYGraph = ROOT.TGraph()
440  extZYGraph.SetMarkerColor(30)
441  extZYGraph.SetMarkerSize(2.25)
442  extZYGraph.SetMarkerStyle(21)
443  extZYSGraph = [0, 0, 0, 0, 0, 0, 0, 0]
444  for sector in range(0, 8):
445  extZYSGraph[sector] = ROOT.TGraph()
446  extZYSGraph[sector].SetMarkerColor(30)
447  extZYSGraph[sector].SetMarkerSize(2.25)
448  extZYSGraph[sector].SetMarkerStyle(21)
449  j = -1
450  for eposition in epositions:
451  j = j + 1
452  x = eposition[0]
453  y = eposition[1]
454  z = eposition[2]
455  extXYGraph.SetPoint(j, x, y)
456  extZYGraph.SetPoint(j, z, y)
457  sector = int(eposition[3])
458  nPoint = extZYSGraph[sector].GetN()
459  extZYSGraph[sector].SetPoint(nPoint, z, abs(x * self.cosine[sector] + y * self.sine[sector]))
460 
461  # Process the MuidHits for event display
462 
463  zMuids = [0, 0, 0, 0, 0, 0, 0, 0]
464  nMuids = [0, 0, 0, 0, 0, 0, 0, 0]
465  muidXYGraph = ROOT.TGraph()
466  muidXYGraph.SetMarkerColor(5)
467  muidXYGraph.SetMarkerSize(2.0)
468  muidXYGraph.SetMarkerStyle(20)
469  muidZYGraph = ROOT.TGraph()
470  muidZYGraph.SetMarkerColor(5)
471  muidZYGraph.SetMarkerSize(2.0)
472  muidZYGraph.SetMarkerStyle(20)
473  muidZYSGraph = [0, 0, 0, 0, 0, 0, 0, 0]
474  for sector in range(0, 8):
475  muidZYSGraph[sector] = ROOT.TGraph()
476  muidZYSGraph[sector].SetMarkerColor(5)
477  muidZYSGraph[sector].SetMarkerSize(2.0)
478  muidZYSGraph[sector].SetMarkerStyle(20)
479  j = -1
480  for muidhit in muidhits:
481  j += 1
482  muidPosition = muidhit.getExtPosition()
483  x = muidPosition[0]
484  y = muidPosition[1]
485  z = muidPosition[2]
486  muidXYGraph.SetPoint(j, x, y)
487  muidZYGraph.SetPoint(j, z, y)
488  sector = muidhit.getSector()
489  nPoint = muidZYSGraph[sector].GetN()
490  muidZYSGraph[sector].SetPoint(nPoint, z, abs(x * self.cosine[sector] + y * self.sine[sector]))
491  if nMuids[sector] == 0:
492  zMuids[sector] = z
493  nMuids[sector] += 1
494 
495  # Process the BKLMHit2ds for event display
496 
497  rpcHitCount = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
498  promptXYGraph = ROOT.TGraph()
499  promptXYGraph.SetMarkerColor(4)
500  promptXYGraph.SetMarkerSize(2.0)
501  promptXYGraph.SetMarkerStyle(29)
502  promptZYGraph = ROOT.TGraph()
503  promptZYGraph.SetMarkerColor(4)
504  promptZYGraph.SetMarkerSize(2.0)
505  promptZYGraph.SetMarkerStyle(29)
506  bkgdXYGraph = ROOT.TGraph()
507  bkgdXYGraph.SetMarkerColor(2)
508  bkgdXYGraph.SetMarkerSize(2.0)
509  bkgdXYGraph.SetMarkerStyle(29)
510  bkgdZYGraph = ROOT.TGraph()
511  bkgdZYGraph.SetMarkerColor(2)
512  bkgdZYGraph.SetMarkerSize(2.0)
513  bkgdZYGraph.SetMarkerStyle(29)
514  promptZYSGraph = [0, 0, 0, 0, 0, 0, 0, 0]
515  bkgdZYSGraph = [0, 0, 0, 0, 0, 0, 0, 0]
516  for sector in range(0, 8):
517  promptZYSGraph[sector] = ROOT.TGraph()
518  promptZYSGraph[sector].SetMarkerColor(4)
519  promptZYSGraph[sector].SetMarkerSize(2.0)
520  promptZYSGraph[sector].SetMarkerStyle(29)
521  bkgdZYSGraph[sector] = ROOT.TGraph()
522  bkgdZYSGraph[sector].SetMarkerColor(2)
523  bkgdZYSGraph[sector].SetMarkerSize(2.0)
524  bkgdZYSGraph[sector].SetMarkerStyle(29)
525  jPrompt = -1
526  jBkgd = -1
527  for hit2d in hit2ds:
528  key = hit2d.getModuleID()
529  layer = (key & self.BKLM_LAYER_MASK) >> self.BKLM_LAYER_BIT
530  sector = (key & self.BKLM_SECTOR_MASK) >> self.BKLM_SECTOR_BIT
531  fb = (key & self.BKLM_SECTION_MASK) >> self.BKLM_SECTION_BIT
532  phiStripMin = hit2d.getPhiStripMin() - 1
533  phiStripMax = hit2d.getPhiStripMax() - 1
534  zStripMin = hit2d.getZStripMin() - 1
535  zStripMax = hit2d.getZStripMax() - 1
536  sectorFB = sector if fb == 0 else sector + 8
537  if layer >= 2:
538  rpcHitCount[sectorFB] += 1
539  dc = self.sectorFBToDC[sectorFB]
540  copper = dc & 0x03
541  finesse = dc >> 2
542  n = rawklms[copper].GetDetectorNwords(0, finesse) >> 1
543  trigCtime = (rawklms[copper].GetTTCtime(0) & 0x07ffffff) << 3
544  ctDiffMax = 99999
545  tCal = -1
546  jZ = -1
547  jPhi = -1
548  ctZ = 0
549  ctPhi = 0
550  for j in range(0, n):
551  if layer != rawLayer[dc][j]:
552  continue
553  if sector != rawSector[dc][j]:
554  continue
555  if fb != rawFb[dc][j]:
556  continue
557  strip = rawStrip[dc][j]
558  plane = rawPlane[dc][j]
559  if plane == 0: # it's a z strip
560  if strip < zStripMin:
561  continue
562  if strip > zStripMax:
563  continue
564  ctZ = rawCtime[dc][j] << 3 # in ns, range is only 8 bits in SCROD (??)
565  jZ = j
566  else: # it's a phi strip
567  if strip < phiStripMin:
568  continue
569  if strip > phiStripMax:
570  continue
571  ctPhi = rawCtime[dc][j] << 3 # in ns, range is only 8 bits in SCROD (??)
572  jPhi = j
573  if (jZ >= 0) and (jPhi >= 0):
574  if layer < 2: # it's a scint layer
575  if abs(ctZ - ctPhi) > 40:
576  continue
577  ct = int((ctZ + ctPhi) * 0.5 - trigCtime - self.ct0Scint[sectorFB]) & 0x3ff
578  if abs(ct - self.ct0Cal) < ctDiffMax:
579  ctDiffMax = int(abs(ct - self.ct0Cal))
580  tCal = ct
581  if ctDiffMax == 0:
582  break
583  else: # it's an RPC layer
584  tCal = ((int(hit2d.getTime()) - trigCtime) & 0x03ff) - self.t0RPC[sectorFB] - 0.75 * jPhi - 0.75 * jZ
585  break
586  x = hit2d.getGlobalPositionX()
587  y = hit2d.getGlobalPositionY()
588  z = hit2d.getGlobalPositionZ()
589  isPromptHit = False
590  if layer < 2:
591  if abs(tCal - self.ct0Cal2d) < 20:
592  isPromptHit = True
593  else:
594  if abs(tCal - self.t0Cal2d) < 20:
595  isPromptHit = True
596  if isPromptHit:
597  jPrompt += 1
598  promptXYGraph.SetPoint(jPrompt, x, y)
599  promptZYGraph.SetPoint(jPrompt, z, y)
600  nPoint = promptZYSGraph[sector].GetN()
601  promptZYSGraph[sector].SetPoint(nPoint, z, abs(x * self.cosine[sector] + y * self.sine[sector]))
602  else:
603  jBkgd += 1
604  bkgdXYGraph.SetPoint(jBkgd, x, y)
605  bkgdZYGraph.SetPoint(jBkgd, z, y)
606  nPoint = bkgdZYSGraph[sector].GetN()
607  bkgdZYSGraph[sector].SetPoint(nPoint, z, abs(x * self.cosine[sector] + y * self.sine[sector]))
608 
609  hasEnoughRPCHits = False
610  for count in rpcHitCount:
611  if count > self.minRPCHits:
612  hasEnoughRPCHits = True
613  break
614  if hasEnoughRPCHits and (len(muidhits) > self.minMuidHits):
615  self.eventDisplays += 1
616  title = 'e{0:02d}r{1}: event {2}'.format(int(self.exp), int(self.run), event)
617  self.hist_XY.SetTitle(title)
618  self.hist_ZY.SetTitle(title)
619  self.eventCanvas.cd(1)
620  self.hist_XY.Draw()
621  for g in self.bklmXY:
622  g.Draw("L")
623  if extXYGraph.GetN() > 0:
624  extXYGraph.Draw("P")
625  if muidXYGraph.GetN() > 0:
626  muidXYGraph.Draw("P")
627  if bkgdXYGraph.GetN() > 0:
628  bkgdXYGraph.Draw("P")
629  if promptXYGraph.GetN() > 0:
630  promptXYGraph.Draw("P")
631  self.eventCanvas.cd(2)
632  self.hist_ZY.Draw()
633  for g in self.bklmZY:
634  g.Draw("L")
635  if extZYGraph.GetN() > 0:
636  extZYGraph.Draw("P")
637  if muidZYGraph.GetN() > 0:
638  muidZYGraph.Draw("P")
639  if bkgdZYGraph.GetN() > 0:
640  bkgdZYGraph.Draw("P")
641  if promptZYGraph.GetN() > 0:
642  promptZYGraph.Draw("P")
643  self.lastTitle = "Title:E{0} (#{1})".format(event, self.eventCounter)
644  self.eventCanvas.Print(self.eventPdfName, self.lastTitle)
645  for sector in range(0, 8):
646  if nMuids[sector] > 0:
647  title = 'e{0:02d}r{1}: event {2} sector {3}'.format(int(self.exp), int(self.run), event, sector)
648  self.hist_XYS[sector].SetTitle(title)
649  self.eventCanvas.cd(1)
650  self.hist_XYS[sector].Draw()
651  for g in self.bklmXY:
652  g.Draw("L")
653  if extXYGraph.GetN() > 0:
654  extXYGraph.Draw("P")
655  if muidXYGraph.GetN() > 0:
656  muidXYGraph.Draw("P")
657  if bkgdXYGraph.GetN() > 0:
658  bkgdXYGraph.Draw("P")
659  if promptXYGraph.GetN() > 0:
660  promptXYGraph.Draw("P")
661  self.eventCanvas.cd(2)
662  z0 = zMuids[sector]
663  self.hist_ZYS.SetTitle(title)
664  self.hist_ZYS.SetBins(10, z0 - 150.0, z0 + 150.0, 10, 125.0, 425.0)
665  self.hist_ZYS.Draw()
666  for g in self.bklmZY:
667  g.Draw("L")
668  for g in self.bklmZYL:
669  g.Draw("L")
670  if extZYSGraph[sector].GetN() > 0:
671  extZYSGraph[sector].Draw("P")
672  if muidZYSGraph[sector].GetN() > 0:
673  muidZYSGraph[sector].Draw("P")
674  if bkgdZYSGraph[sector].GetN() > 0:
675  bkgdZYSGraph[sector].Draw("P")
676  if promptZYSGraph[sector].GetN() > 0:
677  promptZYSGraph[sector].Draw("P")
678  self.lastTitle = "Title:E{0} sector {1}".format(event, sector)
679  self.eventCanvas.Print(self.eventPdfName, self.lastTitle)
EventDisplayer.EventDisplayer.BKLM_SECTION_BIT
int BKLM_SECTION_BIT
bit position for section [0..1]; forward is 0
Definition: EventDisplayer.py:39
EventDisplayer.EventDisplayer.exp
exp
internal copy of experiment number
Definition: EventDisplayer.py:70
EventDisplayer.EventDisplayer.hist_XY
hist_XY
blank scatterplot to define the bounds of the BKLM end view
Definition: EventDisplayer.py:113
EventDisplayer.EventDisplayer.BKLM_ID
int BKLM_ID
COPPER base identifier for BKLM readout.
Definition: EventDisplayer.py:27
EventDisplayer.EventDisplayer.hist_XYS
hist_XYS
list of blank scatterplots to define the per-sector bounds of the BKLM end view
Definition: EventDisplayer.py:125
EventDisplayer.EventDisplayer.lastTitle
lastTitle
title of the last-drawn event display (needed for PDF table of contents' last event)
Definition: EventDisplayer.py:86
EventDisplayer.EventDisplayer.BKLM_SECTOR_BIT
int BKLM_SECTOR_BIT
bit position for sector-1 [0..7]; 0 is on the +x axis and 2 is on the +y axis
Definition: EventDisplayer.py:37
EventDisplayer.EventDisplayer.terminate
def terminate(self)
Definition: EventDisplayer.py:310
EventDisplayer.EventDisplayer.hist_ZY
hist_ZY
blank scatterplot to define the bounds of the BKLM side view
Definition: EventDisplayer.py:116
EventDisplayer.EventDisplayer
Definition: EventDisplayer.py:23
Belle2::PyStoreObj
a (simplified) python wrapper for StoreObjPtr.
Definition: PyStoreObj.h:69
EventDisplayer.EventDisplayer.ct0Cal
ct0Cal
scint-ctime calibration adjustment (ns) for rawKLMs
Definition: EventDisplayer.py:302
EventDisplayer.EventDisplayer.endRun
def endRun(self)
Definition: EventDisplayer.py:322
EventDisplayer.EventDisplayer.t0Cal2d
t0Cal2d
RPC-time calibration adjustment (ns) for BKLMHit2ds.
Definition: EventDisplayer.py:300
EventDisplayer.EventDisplayer.BKLM_PLANE_BIT
int BKLM_PLANE_BIT
bit position for plane-1 [0..1]; 0 is inner-plane
Definition: EventDisplayer.py:33
EventDisplayer.EventDisplayer.eventPdfName
eventPdfName
internal copy of the pathname of the output event-display PDF file
Definition: EventDisplayer.py:74
EventDisplayer.EventDisplayer.run
run
internal copy of run number
Definition: EventDisplayer.py:72
EventDisplayer.EventDisplayer.electIdToModuleId
electIdToModuleId
readout <-> detector map (from the information retrieved from the conditions database)
Definition: EventDisplayer.py:292
EventDisplayer.EventDisplayer.eventDisplays
eventDisplays
event-display counter
Definition: EventDisplayer.py:84
EventDisplayer.EventDisplayer.eventCounter
eventCounter
event counter (needed for PDF table of contents' ordinal event#)
Definition: EventDisplayer.py:82
EventDisplayer.EventDisplayer.bklmXY
bklmXY
list of line-segment (x,y) points for the BKLM end view
Definition: EventDisplayer.py:147
bklmDB.fillDB
def fillDB()
Definition: bklmDB.py:9
EventDisplayer.EventDisplayer.BKLM_SECTOR_MASK
tuple BKLM_SECTOR_MASK
bit mask for sector-1 [0..7]; 0 is on the +x axis and 2 is on the +y axis
Definition: EventDisplayer.py:49
EventDisplayer.EventDisplayer.initialize
def initialize(self)
Definition: EventDisplayer.py:88
EventDisplayer.EventDisplayer.ct0Cal2d
ct0Cal2d
scint-ctime calibration adjustment (ns) for BKLMHit2ds
Definition: EventDisplayer.py:304
EventDisplayer.EventDisplayer.hist_ZYS
hist_ZYS
blank scatterplot to define the per-sector bounds of the rotated BKLM side view
Definition: EventDisplayer.py:143
EventDisplayer.EventDisplayer.sine
sine
table of sines for the BKLM-sector normals
Definition: EventDisplayer.py:107
EventDisplayer.EventDisplayer.ct0Scint
ct0Scint
per-sector variations in scint-ctime calibration adjustment (ns) for rawKLMs
Definition: EventDisplayer.py:308
EventDisplayer.EventDisplayer.beginRun
def beginRun(self)
Definition: EventDisplayer.py:317
EventDisplayer.EventDisplayer.sectorFBToDC
sectorFBToDC
map for sectorFB -> data concentrator
Definition: EventDisplayer.py:294
EventDisplayer.EventDisplayer.BKLM_PLANE_MASK
tuple BKLM_PLANE_MASK
bit mask for plane-1 [0..1]; 0 is inner-plane
Definition: EventDisplayer.py:45
EventDisplayer.EventDisplayer.__init__
def __init__(self, exp, run, eventPdfName, maxDisplays, minRPCHits, minMuidHits)
Definition: EventDisplayer.py:57
EventDisplayer.EventDisplayer.cosine
cosine
table of cosines for the BKLM-sector normals
Definition: EventDisplayer.py:105
EventDisplayer.EventDisplayer.BKLM_LAYER_MASK
tuple BKLM_LAYER_MASK
bit mask for layer-1 [0..15]; 0 is innermost and 14 is outermost
Definition: EventDisplayer.py:47
EventDisplayer.EventDisplayer.eventCanvas
eventCanvas
TCanvas on which event displays will be drawn.
Definition: EventDisplayer.py:96
EventDisplayer.EventDisplayer.BKLM_SECTION_MASK
tuple BKLM_SECTION_MASK
bit mask for section [0..1]; forward is 0
Definition: EventDisplayer.py:51
EventDisplayer.EventDisplayer.event
def event(self)
Definition: EventDisplayer.py:327
EventDisplayer.EventDisplayer.t0RPC
t0RPC
per-sector variations in RPC-time calibration adjustment (ns) for rawKLMs
Definition: EventDisplayer.py:306
Belle2::PyStoreArray
a (simplified) python wrapper for StoreArray.
Definition: PyStoreArray.h:58
EventDisplayer.EventDisplayer.BKLM_STRIP_BIT
int BKLM_STRIP_BIT
bit position for strip-1 [0..47]
Definition: EventDisplayer.py:31
EventDisplayer.EventDisplayer.bklmZYL
bklmZYL
list of line-segment (z,y) points for the BKLM sector's zoomed and rotated side view
Definition: EventDisplayer.py:280
EventDisplayer.EventDisplayer.minMuidHits
minMuidHits
internal copy of the minimum number of MuidHits in the event for event display
Definition: EventDisplayer.py:80
EventDisplayer.EventDisplayer.dcToSectorFB
dcToSectorFB
map for data concentrator -> sectorFB
Definition: EventDisplayer.py:296
EventDisplayer.EventDisplayer.minRPCHits
minRPCHits
internal copy of the minimum number of RPC BKLMHit2ds in any sector for event display
Definition: EventDisplayer.py:78
EventDisplayer.EventDisplayer.BKLM_LAYER_BIT
int BKLM_LAYER_BIT
bit position for layer-1 [0..14]; 0 is innermost
Definition: EventDisplayer.py:35
EventDisplayer.EventDisplayer.t0Cal
t0Cal
RPC-time calibration adjustment (ns) for rawKLMs.
Definition: EventDisplayer.py:298
EventDisplayer.EventDisplayer.bklmZY
bklmZY
list of line-segment (z,y) points for the BKLM side view
Definition: EventDisplayer.py:196
EventDisplayer.EventDisplayer.BKLM_STRIP_MASK
int BKLM_STRIP_MASK
bit mask for strip-1 [0..47]
Definition: EventDisplayer.py:43
EventDisplayer.EventDisplayer.maxDisplays
maxDisplays
internal copy of the maximum number of event displays to write
Definition: EventDisplayer.py:76
EventDisplayer.EventDisplayer.EKLM_ID
int EKLM_ID
COPPER base identifier for EKLM readout.
Definition: EventDisplayer.py:29