Belle II Software  release-05-01-25
test_mixing_overlay_equivalence.py
1 import basf2
2 import ROOT
3 from ROOT import Belle2
4 import b2test_utils
5 from svd import add_svd_simulation
6 
7 # ====================================================================
8 # This test does the following:
9 # 1. Simulate a few hundred ParticleGun events to produce a background
10 # sample.
11 # 2. Produce background overlay sample from background sample. .
12 # 3. Set aside SimHits from a selected sensor to create signal-background
13 # collisions.
14 # 4. Inject collision SimHits, add background using the background mixer,
15 # and digitize. Sort and save.
16 # 5. Inject collision SimHits, digitize, and use background overlay
17 # to add background. Sort and save.
18 # 7. Compare SVDShaperDigits in the two ROOT files.
19 #
20 # NB: I sort both digit samples before saving. Normally, the digits
21 # require sorting only after background overlay, but SVDDigitizer happens
22 # to sort digits differently than ShaperDigitSorter. Both sortings are
23 # correct, but the sorting to eventually correct is that in the
24 # digitizer.
25 # NB: Switching off electronic effects in SVDDigitizer prevents it from
26 # adding noise digits.
27 # ====================================================================
28 
29 
30 xsimhits = []
31 
32 
33 class SetAsideSimHits(basf2.Module):
34  """Set aside all SimHits on a selected sensor/side to create
35  signal/background collisions."""
36 
37  def __init__(self):
38  '''initialize python module'''
39  super().__init__()
40  self.selected_sensorID = Belle2.VxdID(3, 1, 1)
41  '''selected sensor info'''
42 
43  def event(self):
44  '''event'''
45  global xsimhits
46  storesimhits = Belle2.PyStoreArray("SVDSimHits")
47  for h in storesimhits:
48  if h.getSensorID == self.selected_sensorID:
49  xsimhits.append(h)
50 
51 
52 class InjectSimHits(basf2.Module):
53  '''Inject stored SVDSimHits'''
54 
55  def __init__(self):
56  '''initialize python module'''
57  super().__init__()
58  self.simhits = Belle2.PyStoreArray("SVDSimHits")
59  '''sim hit store array'''
60 
61  def initialize(self):
62  '''initialize'''
63  self.simhits.registerInDataStore()
64 
65  def event(self):
66  '''event'''
67 
68  global xsimhits
69  for h in xsimhits:
70  simhit = self.simhits.appendNew()
71  simhit.__assign__(h)
72 
73 
74 if __name__ == "__main__":
75 
77 
78  basf2.B2INFO('Creating background data...')
79 
80  create_bgfile = basf2.create_path()
81  create_bgfile.add_module('EventInfoSetter', expList=[0], runList=[0], evtNumList=[100])
82  create_bgfile.add_module('Gearbox')
83  create_bgfile.add_module('Geometry', components=['MagneticField', 'SVD'])
84  create_bgfile.add_module('ParticleGun')
85  create_bgfile.add_module('FullSim')
86  create_bgfile.add_module('BeamBkgTagSetter', backgroundType="twoPhoton", realTime=310.0)
87  create_bgfile.add_module('RootOutput', outputFileName='bgForMixing.root', branchNames=['BackgroundMetaData', 'SVDSimHits'])
89  b2test_utils.safe_process(create_bgfile)
90 
91  basf2.B2INFO('Creating overlay data...')
92 
93  create_ovrfile = basf2.create_path()
94  create_ovrfile.add_module('EventInfoSetter', expList=[0], runList=[0], evtNumList=[1])
95  create_ovrfile.add_module('Gearbox')
96  create_ovrfile.add_module('Geometry', components=['MagneticField', 'SVD'])
97  create_ovrfile.add_module('BeamBkgMixer', backgroundFiles=['bgForMixing.root'], minTime=-150, maxTime=150)
98  # Turn off generation of noise digits in SVDDigitizer.
99  add_svd_simulation(create_ovrfile)
100  for m in create_ovrfile.modules():
101  if m.name() == "SVDDigitizer":
102  m.param('ElectronicEffects', False)
103 
104  create_ovrfile.add_module(SetAsideSimHits())
105  create_ovrfile.add_module('RootOutput', outputFileName='bgForOverlay.root', branchNames=['SVDShaperDigits'])
107  b2test_utils.safe_process(create_ovrfile)
108 
109  basf2.B2INFO('Background mixing...')
110 
111  produce_mixed = basf2.create_path()
112  produce_mixed.add_module('EventInfoSetter', expList=[0], runList=[0], evtNumList=[1])
113  produce_mixed.add_module('Gearbox')
114  produce_mixed.add_module('Geometry', components=['MagneticField', 'SVD'])
115  # Inject stored simhits
116  produce_mixed.add_module(InjectSimHits())
117  # Time window for background file to just cover one event
118  produce_mixed.add_module('BeamBkgMixer', backgroundFiles=['bgForMixing.root'], minTime=-150, maxTime=150)
119  # Turn off generation of noise digits in SVDDigitizer.
120  add_svd_simulation(produce_mixed)
121  for m in produce_mixed.modules():
122  if m.name() == "SVDDigitizer":
123  m.param('ElectronicEffects', False)
124 
125  produce_mixed.add_module('SVDShaperDigitSorter')
126  produce_mixed.add_module('RootOutput', outputFileName='mixedBg.root')
128  b2test_utils.safe_process(produce_mixed)
129 
130  basf2.B2INFO('Background overlay...')
131 
132  produce_overlaid = basf2.create_path()
133  produce_overlaid.add_module('EventInfoSetter', expList=[0], runList=[0], evtNumList=[1])
134  produce_overlaid.add_module('Gearbox')
135  produce_overlaid.add_module('Geometry', components=['MagneticField', 'SVD'])
136  produce_overlaid.add_module('BGOverlayInput', inputFileNames=['bgForOverlay.root'])
137  produce_overlaid.add_module(InjectSimHits())
138 
139  add_svd_simulation(produce_overlaid)
140  for m in produce_overlaid.modules():
141  if m.name() == "SVDDigitizer":
142  m.param('ElectronicEffects', False)
143 
144  produce_overlaid.add_module('BGOverlayExecutor')
145  # Sort digits after overlay
146  produce_overlaid.add_module('SVDShaperDigitSorter')
147  produce_overlaid.add_module('RootOutput', outputFileName='overlaidBg.root')
149  b2test_utils.safe_process(produce_overlaid)
150 
151  basf2.B2INFO('Comparing...')
152 
153  # Compare diigts stored in the two root files
154  f_mixed = ROOT.TFile('mixedBg.root')
155  tree_mixed = f_mixed.Get('tree')
156 
157  f_over = ROOT.TFile('overlaidBg.root')
158  tree_over = f_over.Get('tree')
159 
160  # We only have a single event in both files
161  n_good = 0
162  for evt1, evt2 in zip(tree_mixed, tree_over):
163  d_mixed = evt1.SVDShaperDigits
164  d_over = evt2.SVDShaperDigits
165  if (len(d_mixed) != len(d_over)):
166  basf2.B2ERROR('Lengths of SVDShaperDigit arrays DIFFER!')
167  break
168  for d1, d2 in zip(d_mixed, d_over):
169  match = True
170  match = match and (d1.getSensorID() == d2.getSensorID())
171  match = match and (d1.isUStrip() == d2.isUStrip())
172  match = match and (d1.getCellID() == d2.getCellID())
173  for s1, s2 in zip(d1.getSamples(), d2.getSamples()):
174  match = match and (int(s1) == int(s2))
175  match = match and (d1.getModeByte() == d2.getModeByte())
176  if not match:
177  print(d1.toString())
178  print(d2.toString())
179  basf2.B2ERROR('Digits do not match.')
180  else:
181  n_good += 1
182 
183  basf2.B2INFO('Processed {0} matching digits.'.format(n_good))
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
test_mixing_overlay_equivalence.SetAsideSimHits.selected_sensorID
selected_sensorID
Definition: test_mixing_overlay_equivalence.py:40
test_mixing_overlay_equivalence.InjectSimHits
Definition: test_mixing_overlay_equivalence.py:52
test_mixing_overlay_equivalence.SetAsideSimHits
Definition: test_mixing_overlay_equivalence.py:33
test_mixing_overlay_equivalence.InjectSimHits.__init__
def __init__(self)
Definition: test_mixing_overlay_equivalence.py:55
b2test_utils.show_only_errors
def show_only_errors()
Definition: __init__.py:87
test_mixing_overlay_equivalence.InjectSimHits.event
def event(self)
Definition: test_mixing_overlay_equivalence.py:65
test_mixing_overlay_equivalence.SetAsideSimHits.__init__
def __init__(self)
Definition: test_mixing_overlay_equivalence.py:37
b2test_utils.clean_working_directory
def clean_working_directory()
Definition: __init__.py:176
Belle2::PyStoreArray
a (simplified) python wrapper for StoreArray.
Definition: PyStoreArray.h:58
b2test_utils.safe_process
def safe_process(*args, **kwargs)
Definition: __init__.py:224
test_mixing_overlay_equivalence.InjectSimHits.initialize
def initialize(self)
Definition: test_mixing_overlay_equivalence.py:61
test_mixing_overlay_equivalence.InjectSimHits.simhits
simhits
Definition: test_mixing_overlay_equivalence.py:58
test_mixing_overlay_equivalence.SetAsideSimHits.event
def event(self)
Definition: test_mixing_overlay_equivalence.py:43