Belle II Software  release-05-01-25
__init__.py
1 import basf2 as b2
2 from basf2.pickle_path import deserialize_path, deserialize_module, serialize_path, serialize_module
3 
4 import sys
5 import ROOT
6 from ROOT import Belle2
7 from caf.framework import Calibration
8 
9 from collections import OrderedDict
10 
11 b2.set_log_level(b2.LogLevel.INFO)
12 UVWABC = [1, 2, 3, 4, 5, 6]
13 
14 
15 def collect(calibration, basf2_args=None):
16  import pickle
17  import subprocess
18 
19  if basf2_args is None:
20  basf2_args = []
21 
22  main = b2.create_path()
23  main.add_module('RootInput', inputFileNames=calibration.input_files)
24  main.add_module('HistoManager', histoFileName='RootOutput.root')
25  main.add_path(calibration.pre_collector_path)
26  main.add_module(calibration.collector)
27 
28  path_file_name = calibration.name + '.path'
29  with open(path_file_name, 'bw') as serialized_path_file:
30  pickle.dump(serialize_path(main), serialized_path_file)
31 
32  subprocess.call(["basf2", "--execute-path", path_file_name] + basf2_args)
33 
34 
35 def calibrate(calibration, input_file='RootOutput.root'):
36  for algo in calibration.algorithms:
37  algo.algorithm.setInputFileNames([input_file])
38  algo.algorithm.execute()
39  algo.algorithm.commit()
40 
41 
43  """ The generic Millepede calibration collector+algorithm """
44 
45  def __init__(self,
46  components=None,
47  tracks=None,
48  particles=None,
49  vertices=None,
50  primary_vertices=None,
51  path=None):
52  """
53  components are the names of DB objects to calibrate (BeamSpot etc.)
54  tracks are collections of RecoTracks of fitted charged tracks (usually cosmic rays, no associated particle)
55  particles are names of ParticleLists with single charged particles
56  vertices are names ParticleLists with at least two body decays fitted with vertex constraint
57  primary_vertices are names of ParticleLists with at least two body decays fitted with beam+vertex constraint
58  """
59 
60  if components is None:
61  components = []
62 
63  # By default, at least the raw RecoTracks are fitted (let's specify its name explicitly here)
64  if (tracks is None) and (particles is None) and (vertices is None) and (primary_vertices is None):
65  tracks = ['RecoTracks']
66 
67  if particles is None:
68  particles = []
69 
70  if vertices is None:
71  vertices = []
72 
73  if primary_vertices is None:
74  primary_vertices = []
75 
76  if path is None:
77  path = b2.create_path()
78  path.add_module('Gearbox')
79  path.add_module('Geometry')
80  path.add_module('SetupGenfitExtrapolation', noiseBetheBloch=False, noiseCoulomb=False, noiseBrems=False)
81 
82 
84 
85  self.path = path
86 
87  self.collector = b2.register_module('MillepedeCollector')
88 
89  self.parameters = OrderedDict()
90 
91  self.commands = OrderedDict()
92 
93  self.set_components(components)
94 
95  self.collector.param('tracks', tracks)
96  self.collector.param('particles', particles)
97  self.collector.param('vertices', vertices)
98  self.collector.param('primaryVertices', primary_vertices)
99  self.collector.param('minPValue', 0.)
100  self.collector.param('useGblTree', True)
101 
102  self.set_command('method diagonalization 3 0.1')
103  self.set_command('skipemptycons')
104  self.set_command('entries 10')
105  self.set_command('hugecut 50')
106  self.set_command('chiscut 30. 6.')
107  self.set_command('outlierdownweighting 3')
108  self.set_command('dwfractioncut 0.1')
109 
110  def set_components(self, components):
111  """ Select db objects for calibration from list of their names """
112 
113  self.components = components
114  calibrate_vertex = ((components == []) or ('BeamSpot' in components))
115 
116  std_components = ROOT.vector('string')()
117  for component in components:
118  std_components.push_back(component)
119 
120  self.algo.setComponents(std_components)
121  self.collector.param('components', components)
122  self.collector.param('calibrateVertex', calibrate_vertex)
123 
124  def set_command(self, cmd_name, command=''):
125  """ Set command for Millepede steering """
126  if command is None:
127  self.commands[cmd_name] = ''
128  return
129 
130  if command == '':
131  command = cmd_name
132  else:
133  command = cmd_name + ' ' + str(command)
134 
135  cmd_words = str(command).split(' ')
136  self.commands[cmd_words[0]] = command
137 
138  def pre_algo(self, algorithm, iteration):
139  """ Fcn to execute before algorithm... """
140  pass
141 
142  def create(self, name, input_files):
143  """ Create the CAF Calibration object """
144 
145  # We keep steering options for pede in python
146  # as long as possible for the user before passing to algo
147  for command in self.commands.values():
148  self.algo.steering().command(command)
149  # Add parameter fixing/setting
150  self.algo.steering().command('Parameters')
151  for command in self.parameters.values():
152  self.algo.steering().command(command)
153 
154  # If you call this, you are going to use calibration framework most likely.
155  # This NOW does not anymore move/copy files from collector to algorithm, so
156  # we have to remeber where they are located at the time of creation
157  if not self.get_param('useGblTree'):
158  self.set_param(True, 'absFilePaths')
159 
160  cal = Calibration(name, self.collector, self.algo, input_files)
161  cal.output_patterns.append('Belle2FileCatalog.xml')
162  cal.output_patterns.append('constraints.txt')
163  cal.output_patterns.append('*.mille')
164  # Adding in setup basf2 paths and functions for the collector and algorithm
165  cal.pre_collector_path = self.path
166  cal.pre_algorithms = self.pre_algo
167 
168  return cal
169 
170  def get_param(self, param, module=None):
171  """ Get parameter of the collector module or any module in path (given its name) """
172  module = self.get_module(module)
173  for mpi in module.available_params():
174  if mpi.name == param:
175  return mpi.values
176  return None
177 
178  def set_param(self, value, param, module=None):
179  """ Set parameter of the collector module or any module in path (given its name) """
180  module = self.get_module(module)
181  if module is not None:
182  module.param(param, value)
183 
184  def get_module(self, module=None):
185  """ Get collector module or any other module from path (given its name) """
186  if module is None or module == self.collector.name() or module == self.collector:
187  module = self.collector
188  else:
189  for mod in self.path.modules():
190  if mod.name() == module:
191  module = mod
192  break
193  if isinstance(module, b2.Module):
194  return module
195  return None
196 
197  def __getstate__(self):
198  """ serialization """
199 
200  state = self.__dict__.copy()
201  state['path'] = serialize_path(self.path)
202  state['collector'] = serialize_module(self.collector)
203  return state
204 
205  def __setstate__(self, state):
206  """ de-serialization """
207 
208  self.__dict__.update(state)
209  self.collector = deserialize_module(state['collector'])
210  self.path = deserialize_path(state['path'])
211 
212  def fixGlobalParam(self, uid, element, param, value=0.):
213  """ Generic fcn to manipulate (fix to given value) parameters identified by db object id (uid),
214  element and parameter number """
215  label = Belle2.GlobalLabel()
216  label.construct(uid, element, param)
217  self.parameters[str(label.label())] = (str(label.label()) + ' ' + str(value) + ' -1.')
218 
219  # CDC helpers --------------------------------------------------------------------------------------------------
220 
221  def fixCDCLayerX(self, layer):
222  """ fix CDC layer X-shift """
223  self.fixGlobalParam(Belle2.CDCLayerAlignment.getGlobalUniqueID(), layer, Belle2.CDCLayerAlignment.layerX)
224 
225  def fixCDCLayerY(self, layer):
226  """ fix CDC layer Y-shift """
227  self.fixGlobalParam(Belle2.CDCLayerAlignment.getGlobalUniqueID(), layer, Belle2.CDCLayerAlignment.layerY)
228 
229  def fixCDCLayerRot(self, layer):
230  """ fix CDC layer Phi rotation """
231  self.fixGlobalParam(Belle2.CDCLayerAlignment.getGlobalUniqueID(), layer, Belle2.CDCLayerAlignment.layerPhi)
232 
233  def fixCDCTimeWalk(self, board_id):
234  """ fix CDC time walk parameter for given board """
236 
237  def fixCDCTimeZero(self, wire_id):
238  """ fix CDC tie zero parameter for given wire """
240 
241  # VXD helpers --------------------------------------------------------------------------------------------------
242 
243  def fixVXDid(self, layer, ladder, sensor, segment=0, parameters=[1, 2, 3, 4, 5, 6]):
244  """ Fix VXD element parameters by layer,ladder,sensor and segment number """
245  for param in parameters:
247  layer, ladder, sensor, segment).getID(), param)
248  """
249  {{"PXD.Ying"}, {Belle2::VxdID(1, 0, 0, 1)}},
250  {{"PXD.Yang"}, {Belle2::VxdID(1, 0, 0, 2)}},
251  {{"SVD.Pat"}, {Belle2::VxdID(3, 0, 0, 1)}},
252  {{"SVD.Mat"}, {Belle2::VxdID(3, 0, 0, 2)}}
253  """
254 
255  def fixPXDYing(self, parameters=UVWABC):
256  """ fix PXD Ying half-shell """
257  self.fixVXDid(1, 0, 0, 1, parameters)
258 
259  def fixPXDYang(self, parameters=UVWABC):
260  """ fix PXD Yang half-shell """
261  self.fixVXDid(1, 0, 0, 2, parameters)
262 
263  def fixSVDPat(self, parameters=UVWABC):
264  """ fix SVD Pat half-shell """
265  self.fixVXDid(3, 0, 0, 1, parameters)
266 
267  def fixSVDMat(self, parameters=UVWABC):
268  """ fix SVD Mat half-shell """
269  self.fixVXDid(3, 0, 0, 2, parameters)
270 
271  # EKLM helpers --------------------------------------------------------------------------------------------------
272 
273  def fixEKLMSector(self, endcap, layer, sector, parameters=[1, 2, 6]):
274  """ Fix EKLM sector parameters """
275 
276  for ipar in parameters:
278  Belle2.EKLMElementID(endcap, layer, sector).getGlobalNumber(),
279  ipar)
280 
281  def fixEKLMSegment(self, endcap, layer, sector, plane, segment, parameters=[1, 2, 6]):
282  """ Fix EKLM segment parameters """
283  for ipar in parameters:
285  Belle2.EKLMElementID(endcap, layer, sector, plane, segment).getGlobalNumber(),
286  ipar)
287 
288  # BKLM helpers --------------------------------------------------------------------------------------------------
289 
290  def fixBKLMModule(self, sector, layer, forward, parameters=UVWABC):
291  """ Fix a BKLM module """
292  for ipar in parameters:
293  self.fixGlobalParam(
295  Belle2.BKLMElementID(
296  forward,
297  sector,
298  layer).getID(),
299  ipar)
300 
301  def fixBKLM(self, sectors=range(1, 9), layers=range(1, 16), forbackwards=[0, 1]):
302  """ Fix (all by default) BKLM modules """
303  for sector in sectors:
304  for layer in layers:
305  for forward in forbackwards:
306  self.fixBKLMModule(forward, sector, layer)
alignment.MillepedeCalibration.parameters
parameters
Parameter config at algo level (fixing)
Definition: __init__.py:83
alignment.MillepedeCalibration.fixCDCLayerY
def fixCDCLayerY(self, layer)
Definition: __init__.py:225
Belle2::MillepedeAlgorithm
Class implementing Millepede calibration algorithm.
Definition: MillepedeAlgorithm.h:35
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
alignment.MillepedeCalibration.algo
algo
the algorithm
Definition: __init__.py:77
alignment.MillepedeCalibration.__setstate__
def __setstate__(self, state)
Definition: __init__.py:205
modules
Definition: modules.py:1
alignment.MillepedeCalibration.fixCDCLayerX
def fixCDCLayerX(self, layer)
Definition: __init__.py:221
Belle2::getID
int getID(const std::vector< double > &breaks, double t)
get id of the time point t
Definition: calibTools.h:71
Belle2::VXDAlignment::getGlobalUniqueID
static unsigned short getGlobalUniqueID()
Get global unique id.
Definition: VXDAlignment.h:57
Belle2::BKLMAlignment::getGlobalUniqueID
static unsigned short getGlobalUniqueID()
Get global unique identifier.
Definition: BKLMAlignment.h:73
Belle2::CDCLayerAlignment::getGlobalUniqueID
static unsigned short getGlobalUniqueID()
Get global unique id.
Definition: CDCLayerAlignment.h:58
basf2.pickle_path
Definition: pickle_path.py:1
Belle2::GlobalLabel
Class to convert to/from global labels for Millepede II to/from detector & parameter identificators.
Definition: GlobalLabel.h:51
alignment.MillepedeCalibration.fixGlobalParam
def fixGlobalParam(self, uid, element, param, value=0.)
Definition: __init__.py:212
alignment.MillepedeCalibration.pre_algo
def pre_algo(self, algorithm, iteration)
Definition: __init__.py:138
alignment.MillepedeCalibration.set_command
def set_command(self, cmd_name, command='')
Definition: __init__.py:124
alignment.MillepedeCalibration.__init__
def __init__(self, components=None, tracks=None, particles=None, vertices=None, primary_vertices=None, path=None)
Definition: __init__.py:45
alignment.MillepedeCalibration.fixPXDYing
def fixPXDYing(self, parameters=UVWABC)
Definition: __init__.py:255
alignment.MillepedeCalibration.components
components
db objects for calibration
Definition: __init__.py:113
alignment.MillepedeCalibration.path
path
pre-collector path
Definition: __init__.py:79
alignment.MillepedeCalibration.fixCDCLayerRot
def fixCDCLayerRot(self, layer)
Definition: __init__.py:229
alignment.MillepedeCalibration.get_param
def get_param(self, param, module=None)
Definition: __init__.py:170
alignment.MillepedeCalibration.fixBKLMModule
def fixBKLMModule(self, sector, layer, forward, parameters=UVWABC)
Definition: __init__.py:290
alignment.MillepedeCalibration.fixPXDYang
def fixPXDYang(self, parameters=UVWABC)
Definition: __init__.py:259
alignment.MillepedeCalibration.fixVXDid
def fixVXDid(self, layer, ladder, sensor, segment=0, parameters=[1, 2, 3, 4, 5, 6])
Definition: __init__.py:243
Belle2::CDCTimeZeros::getGlobalUniqueID
static unsigned short getGlobalUniqueID()
Get global unique id.
Definition: CDCTimeZeros.h:150
alignment.MillepedeCalibration.get_module
def get_module(self, module=None)
Definition: __init__.py:184
alignment.MillepedeCalibration.set_components
def set_components(self, components)
Definition: __init__.py:110
alignment.MillepedeCalibration
Definition: __init__.py:42
alignment.MillepedeCalibration.fixEKLMSegment
def fixEKLMSegment(self, endcap, layer, sector, plane, segment, parameters=[1, 2, 6])
Definition: __init__.py:281
Belle2::EKLMAlignment::getGlobalUniqueID
static unsigned short getGlobalUniqueID()
Get global unique identifier.
Definition: EKLMAlignment.h:73
alignment.MillepedeCalibration.collector
collector
the collector module
Definition: __init__.py:81
alignment.MillepedeCalibration.fixSVDPat
def fixSVDPat(self, parameters=UVWABC)
Definition: __init__.py:263
alignment.MillepedeCalibration.commands
commands
Commands for pede steering.
Definition: __init__.py:85
alignment.MillepedeCalibration.__getstate__
def __getstate__(self)
Definition: __init__.py:197
steering
Definition: steering.py:1
alignment.MillepedeCalibration.fixEKLMSector
def fixEKLMSector(self, endcap, layer, sector, parameters=[1, 2, 6])
Definition: __init__.py:273
alignment.MillepedeCalibration.fixCDCTimeWalk
def fixCDCTimeWalk(self, board_id)
Definition: __init__.py:233
alignment.MillepedeCalibration.fixSVDMat
def fixSVDMat(self, parameters=UVWABC)
Definition: __init__.py:267
alignment.MillepedeCalibration.set_param
def set_param(self, value, param, module=None)
Definition: __init__.py:178
Calibration
Definition: Calibration.py:1
alignment.MillepedeCalibration.create
def create(self, name, input_files)
Definition: __init__.py:142
alignment.MillepedeCalibration.fixBKLM
def fixBKLM(self, sectors=range(1, 9), layers=range(1, 16), forbackwards=[0, 1])
Definition: __init__.py:301
alignment.MillepedeCalibration.fixCDCTimeZero
def fixCDCTimeZero(self, wire_id)
Definition: __init__.py:237
Belle2::CDCTimeWalks::getGlobalUniqueID
static unsigned short getGlobalUniqueID()
Get global unique id.
Definition: CDCTimeWalks.h:165