Belle II Software development
constraints.py
1
8
9import basf2
10from ROOT import Belle2
11
12import os
13import math
14
15
16class Constraint():
17 """
18 Class representing a linear constraint for global parameters
19 """
20
21 def __init__(self, comment="", value=0.):
22 """
23 Initialize constraint
24
25 comment : str
26 TODO: Not yet used
27 value : float
28 The constant term of the constraint: sum(c_i * par_i) = value
29 """
30
31 self.value = value
32
33 self.comment = comment
34
35 self.data = []
36
37 def add(self, label, coeff):
38 """
39 Add coefficient for a global parameter
40
41 Parameters
42 ----------
43 label : int
44 global parameter id (GlobalLabel::label())
45 coeff : float
46 coefficient of the parameter in constraint
47 """
48 entry = (label, coeff)
49 self.data.append(entry)
50
51 def get_checksum(self):
52 """Get a checksum to distinguish constraints quickly (computed only from labels)
53 """
54 labels_only = [label for (label, coeff) in self.data]
55 checksum = hash(str(labels_only))
56 return checksum
57
58
60 """
61 Base class representing a "generator" for file with set of constraints
62
63 Can be used directly if you already have a file with constraints.
64 """
65
66 def __init__(self, filename):
67 """
68 filename : str
69 Desired filename for the constraint file
70 """
71
72 self.filename = filename
73
74 def generate(self):
75 """
76 Should be overridden by the deriving classes and actually
77 fill the dictionary with the constraints (this runs within a basf2 module
78 event function - so you have geometry available)
79 """
80 consts = []
81 return consts
82
83 def configure_collector(self, collector):
84 """
85 Can be overridden be child classes to pass additional configuration to the
86 MillepedeCollector (activated by the use of the constraints)
87 """
88 pass
89
90
91def generate_constraints(constraint_sets, timedep, global_tags, init_event):
92 """
93 Run the constraints generation and return file names
94
95 Parameters
96 ----------
97 constraint_sets : list (alignment.Constraints)
98 List of sets of constraints
99 timedep_config : list(tuple(list(int), list(tuple(int, int, int))))
100 Time-dependence configuration.
101 Each list item is 2-tuple with list of parameter numbers (use alignment.parameters to get them) and
102 the (event, run, exp) numbers at which values of these parameters can change.
103 global_tags : list (str)
104 List of global tag names and/or (absolute) file paths to local databases
105 init_event : tuple( int, int, int)
106 Event (event, run, exp) at which to initialize time-INdependent constraints
107 """
108 files = []
109 for filename in [consts.filename for consts in constraint_sets]:
110 files.append(os.path.abspath(filename))
111
112 from alignment.constraints_generator import save_config
113 ccfn = save_config(constraint_sets, timedep, global_tags, init_event)
114 os.system('basf2 {} {}'.format(Belle2.FileSystem.findFile('alignment/scripts/alignment/constraints_generator.py'), ccfn))
115
116 return files
117
118# ----------------------------- Sub-detector specific constraint classes -----------------------------------------------
119
120
122 """
123 Constraints for VXD hierarchy
124
125 They come in 3 types:
126
127 0 : no hierarchy, just sensors aligned in their local system
128 1 : sensors -> half-shells (no ladders)
129 2 : sensors -> ladders ->half-shells (aka full hierarchy)
130
131 and additionally can use them for PXD, SVD or both (default) (setting pxd/svd to False
132 will disable the hierarchy for that sub-detector)
133
134 NOTE: the filename cannot be currently changed, it is hard-coded in the collector module
135 as these constraints are actually generated by C++ code and not by python modules as others.
136 """
137
138 def __init__(self, type=2, pxd=True, svd=True):
139 """ Initialize
140 type : int
141 Type of constraints (0 - no, 1 - half-shells, 2 -full)
142 pxd : bool
143 Use constraints for PXD?
144 svd : bool
145 Use constraints for SVD?
146 """
147
148 # TODO: cannot currently change the filename in collector, so fixed here
149 super().__init__("constraints.txt")
150
151 self.type = type
152
153 self.pxd = pxd
154
155 self.svd = svd
156
157 def configure_collector(self, collector):
158 """
159 Propagate the hierarchy configuration to the collector
160 """
161 collector.param('hierarchyType', self.type)
162 collector.param('enablePXDHierarchy', self.pxd)
163 collector.param('enableSVDHierarchy', self.svd)
164
165 def generate(self):
166 """
167 Generate the constraints - does nothing. The collector will make the job
168 """
169 print("Generating constraints for VXD (no-op, file created by collector) ...")
170 return []
171
172
174 """
175 Various constraints for CDC layers
176
177 (Code from Claus Kleinwort)
178
179 """
180
181 cdc = [['A1', 8, 160, 16.80, 23.80, 0., -35.9, 67.9],
182 ['U2', 6, 160, 25.70, 34.80, 0.068, -51.4, 98.6],
183 ['A3', 6, 192, 36.52, 45.57, 0., -57.5, 132.9],
184 ['V4', 6, 224, 47.69, 56.69, -0.060, -59.6, 144.7],
185 ['A5', 6, 256, 58.41, 67.41, 0., -61.8, 146.8],
186 ['U6', 6, 288, 69.53, 78.53, 0.064, -63.9, 148.9],
187 ['A7', 6, 320, 80.25, 89.25, 0., -66.0, 151.0],
188 ['V8', 6, 352, 91.37, 100.37, -0.072, -68.2, 153.2],
189 ['A9', 6, 384, 102.00, 111.14, 0., -70.2, 155.3]]
190
191
192 parameter = [(1, 'x-offset bwd'), (2, 'y-offset bwd'), (6, 'z-rotation bwd'),
193 (11, 'x-offset fwd-bwd'), (12, 'y-offset fwd-bwd'), (16, 'z-rotation fwd-bwd')]
194
195 def __init__(self, filename='cdc-constraints.txt', rigid=True, twist=True, z_offset=False, r_scale=False, z_scale=False):
196 """
197 filename : str
198 Filename with constraints (probably generated later)
199 rigid : bool
200 5D CDC constraints: x+y+rot + dx+dy
201 twist : bool
202 CDC twist constraint: drot
203 z_offset : bool
204 Constraint for Z-offset
205 r_scale : bool
206 Constraint for R-scale
207 z_scale : bool
208 Constraint for Z-scale
209 """
210 super().__init__(filename)
211
212 self.rigid = rigid
213
214 self.twist = twist
215
216 self.z_offset = z_offset
217
218 self.r_scale = r_scale
219
220 self.z_scale = z_scale
221
222 def generate(self):
223 """Generate constraints from CDC geometry
224 """
225 print("Generating constraints for CDC layers...")
226
227 def cdc_layer_label(layer, param):
228 wire = 511
229 wireid = Belle2.WireID(layer, wire).getEWire()
230 label = Belle2.GlobalLabel()
231 label.construct(Belle2.CDCAlignment.getGlobalUniqueID(), wireid, param)
232 return label.label()
233
234 def cmp(a, b):
235 return (a > b) - (a < b)
236
237 consts = []
238
239 for par in self.parameter:
240 nTot = 0
241 for sl in self.cdc:
242 nlyr = sl[1]
243 for lyr in range(nlyr):
244 nTot += nlyr
245 nLayer = nTot
246
247 if self.rigid:
248 # all layers
249 for par in self.parameter:
250 if par[0] == 16:
251 continue
252 const = Constraint()
253
254 nTot = 0
255 for sl in self.cdc:
256 nlyr = sl[1]
257 for lyr in range(nlyr):
258 label = cdc_layer_label(nTot + lyr, par[0])
259 const.add(label, 1.0)
260 nTot += nlyr
261
262 consts.append(const)
263
264 if self.twist:
265 # all layers
266 for par in self.parameter:
267 if par[0] != 16:
268 continue
269 const = Constraint()
270
271 nTot = 0
272 for sl in self.cdc:
273 nlyr = sl[1]
274 for lyr in range(nlyr):
275 label = cdc_layer_label(nTot + lyr, par[0])
276 const.add(label, 1.0)
277 nTot += nlyr
278
279 # f.write('\n')
280 consts.append(const)
281
282 if self.z_offset:
283 # stereo layers (Z offset)
284 par = self.parameter[2]
285 const = Constraint()
286
287 nTot = 0
288 for sl in self.cdc:
289 nlyr = sl[1]
290 rInner = sl[3]
291 rOuter = sl[4]
292 stereo = sl[5]
293 if stereo != 0.:
294 dr = (rOuter - rInner) / (nlyr - 1)
295 rWire = rInner
296 for lyr in range(nlyr):
297 label = cdc_layer_label(nTot + lyr, par[0])
298 const.add(label, stereo * rWire)
299 rWire += dr
300 nTot += nlyr
301
302 # f.write('\n')
303 consts.append(const)
304
305 if self.r_scale:
306 # all layers 2nd
307 for par in self.parameter[:2]:
308 const = Constraint()
309
310 nTot = 0
311 for sl in self.cdc:
312 nlyr = sl[1]
313 for lyr in range(nlyr):
314 label = cdc_layer_label(nTot + lyr, par[0])
315 der = cmp(2. * float(nTot + lyr) + 0.5, float(nLayer))
316 const.add(label, der)
317 nTot += nlyr
318
319 consts.append(const)
320
321 if self.z_scale:
322 # stereo layers (Z -scale)
323 par = self.parameter[5]
324 const = Constraint()
325 nTot = 0
326 for sl in self.cdc:
327 nlyr = sl[1]
328 stereo = sl[5]
329 if stereo != 0.:
330 for lyr in range(nlyr):
331 label = cdc_layer_label(nTot + lyr, par[0])
332 const.add(label, stereo)
333 nTot += nlyr
334
335 # f.write('\n')
336 consts.append(const)
337
338 return consts
339
340
342 """
343 Constraint to fix sum of corrections to CDC T0's (per wire) to zero
344 """
345
346
347 wires_in_layer = [
348 160, 160, 160, 160, 160, 160, 160, 160,
349 160, 160, 160, 160, 160, 160,
350 192, 192, 192, 192, 192, 192,
351 224, 224, 224, 224, 224, 224,
352 256, 256, 256, 256, 256, 256,
353 288, 288, 288, 288, 288, 288,
354 320, 320, 320, 320, 320, 320,
355 352, 352, 352, 352, 352, 352,
356 384, 384, 384, 384, 384, 384]
357
358 def __init__(self, filename='cdc-T0-constraints.txt'):
359 """ Initialize
360 filename : str
361 Can use different filename
362 """
363 super().__init__(filename)
364 pass
365
366 def generate(self):
367 """
368 Generate the constraints
369 """
370 print("Generating constraints for CDC T0's ...")
371 consts = []
372 const = Constraint()
373
374 for layer in range(0, 56):
375 for wire in range(0, self.wires_in_layer[layer]):
376 label = Belle2.GlobalLabel()
377 label.construct(Belle2.CDCTimeZeros.getGlobalUniqueID(), Belle2.WireID(layer, wire).getEWire(), 0)
378 const.add(label.label(), 1.0)
379
380 consts.append(const)
381 return consts
382
383
385 """
386 Various constraints for CDC wires (w.r.t. layers)
387 """
388
389
390 wires_in_layer = [
391 160, 160, 160, 160, 160, 160, 160, 160,
392 160, 160, 160, 160, 160, 160,
393 192, 192, 192, 192, 192, 192,
394 224, 224, 224, 224, 224, 224,
395 256, 256, 256, 256, 256, 256,
396 288, 288, 288, 288, 288, 288,
397 320, 320, 320, 320, 320, 320,
398 352, 352, 352, 352, 352, 352,
399 384, 384, 384, 384, 384, 384]
400
402 self,
403 filename='cdc-wire-constraints.txt',
404 layers=None,
405 layer_rigid=True,
406 layer_radius=None,
407 cdc_radius=False,
408 hemisphere=None):
409 """Initialize constraint
410
411 Parameters
412 ----------
413 filename : str
414 Can override the default file name
415 layers : list (int)
416 List of layer numbers for which to generate the constraints (default is 0..55)
417 layer_rigid : bool
418 6 constraints - fix sum of shifts of wire ends in x/y/rotation at fwd/bwd end-plate
419 layer_radius : list(int)
420 2 constraints to fix average change of layer radius (wires in layer moving away from origin) for all layers in list
421 cdc_radius : bool
422 1 constraint - fix average change in CDC radius from all wires
423 hemisphere : list(int)
424 Modifies rigid and layer_radius (layer_radius constraint is added if not present for selected layer(s))
425 constraint by splitting the constraints for a given list of layers to L/R up/down half of CDC
426
427 """
428
429 super().__init__(filename)
430
431 if layers is None:
432 layers = list(range(0, 56))
433
434 self.layers = layers
435
439 self.layer_rigid = layer_rigid
440
442 if layer_radius is None:
443 layer_radius = []
444
445 self.layer_radius = layer_radius
446
448 self.cdc_radius = cdc_radius
449
450 if hemisphere is None:
451 hemisphere = []
452
453 self.hemisphere = hemisphere
454
455 def configure_collector(self, collector):
456 """Enables wire-by-wire derivatives in collector
457 """
458 basf2.B2WARNING("Adding CDC wire constraints -> enabling wire-by-wire alignment derivatives")
459 collector.param('enableWireByWireAlignment', True)
460
461 def get_label(self, layer, wire, parameter):
462 """Return GlobalLabel for wire parameter"""
463
464 wireid = Belle2.WireID(layer, wire).getEWire()
465 label = Belle2.GlobalLabel()
466 label.construct(Belle2.CDCAlignment.getGlobalUniqueID(), wireid, parameter)
467 return label.label()
468
469 def generate(self):
470 """Generate constraints """
471 print("Generating constraints for CDC wires...")
472 consts = []
473
474 layers = self.layers
475
476 if self.layer_rigid:
477 for layer in layers:
478 const = Constraint()
479 # sum of wire X (BWD) in layer
480 for wire in range(0, self.wires_in_layer[layer]):
481 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdX), 1.)
482 consts.append(const)
483
484 for layer in layers:
485 const = Constraint()
486 # sum of wire Y (BWD) in layer
487 for wire in range(0, self.wires_in_layer[layer]):
488 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdY), 1.)
489 consts.append(const)
490
491 for layer in set(layers) - set(self.hemisphere):
492 const = Constraint()
493 # sum of wire rotations (BWD) in layer
494 for wire in range(0, self.wires_in_layer[layer]):
495
496 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getBackwardPos3D().phi()
497
498 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdX), -math.sin(wirePhi))
499 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdY), +math.cos(wirePhi))
500 consts.append(const)
501
502 for layer in layers:
503 const = Constraint()
504 # sum of wire X (FWD) in layer
505 for wire in range(0, self.wires_in_layer[layer]):
506 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdX), 1.)
507 consts.append(const)
508
509 for layer in layers:
510 const = Constraint()
511 # sum of wire Y (FWD) in layer
512 for wire in range(0, self.wires_in_layer[layer]):
513 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdY), 1.)
514 consts.append(const)
515
516 for layer in set(layers) - set(self.hemisphere):
517 const = Constraint()
518 # sum of wire rotations (FWD) in layer
519 for wire in range(0, self.wires_in_layer[layer]):
520
521 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getForwardPos3D().phi()
522
523 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdX), -math.sin(wirePhi))
524 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdY), +math.cos(wirePhi))
525 consts.append(const)
526
527 for layer in self.hemisphere:
528 const = Constraint()
529 # sum of wire rotations (BWD) in layer ... RIGHT hemisphere
530 for wire in range(0, self.wires_in_layer[layer]):
531
532 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getBackwardPos3D().phi()
533
534 if math.cos(wirePhi) <= 0.:
535 continue
536
537 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdX), -math.sin(wirePhi))
538 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdY), +math.cos(wirePhi))
539 consts.append(const)
540
541 for layer in self.hemisphere:
542 const = Constraint()
543 # sum of wire rotations (BWD) in layer ... LEFT hemisphere
544 for wire in range(0, self.wires_in_layer[layer]):
545
546 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getBackwardPos3D().phi()
547
548 if math.cos(wirePhi) > 0.:
549 continue
550
551 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdX), -math.sin(wirePhi))
552 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdY), +math.cos(wirePhi))
553 consts.append(const)
554
555 for layer in self.hemisphere:
556 const = Constraint()
557 # sum of wire rotations (FWD) in layer ... RIGHT
558 for wire in range(0, self.wires_in_layer[layer]):
559
560 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getForwardPos3D().phi()
561
562 if math.cos(wirePhi) <= 0.:
563 continue
564
565 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdX), -math.sin(wirePhi))
566 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdY), +math.cos(wirePhi))
567 consts.append(const)
568
569 for layer in self.hemisphere:
570 const = Constraint()
571 # sum of wire rotations (FWD) in layer ... LEFT
572 for wire in range(0, self.wires_in_layer[layer]):
573
574 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getForwardPos3D().phi()
575
576 if math.cos(wirePhi) > 0.:
577 continue
578
579 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdX), -math.sin(wirePhi))
580 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdY), +math.cos(wirePhi))
581 consts.append(const)
582
583 # layer_radius:
584 for layer in [lyr for lyr in self.layer_radius if lyr not in self.hemisphere]:
585 const = Constraint()
586 # sum of wire rotations (BWD) in layer
587 for wire in range(0, self.wires_in_layer[layer]):
588
589 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getBackwardPos3D().phi()
590
591 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdX), +math.cos(wirePhi))
592 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdY), +math.sin(wirePhi))
593 const.add(self.get_label(layer, 0, 0), 0.)
594 consts.append(const)
595
596 for layer in [lyr for lyr in self.layer_radius if lyr not in self.hemisphere]:
597 const = Constraint()
598 # sum of wire rotations (FWD) in layer
599 for wire in range(0, self.wires_in_layer[layer]):
600
601 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getForwardPos3D().phi()
602
603 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdX), +math.cos(wirePhi))
604 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdY), +math.sin(wirePhi))
605 const.add(self.get_label(layer, 0, 0), 0.)
606 consts.append(const)
607
608 # hemisphere:
609 for layer in self.hemisphere:
610 const = Constraint()
611 # sum of wire rotations (BWD) in layer
612 for wire in range(0, self.wires_in_layer[layer]):
613
614 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getBackwardPos3D().phi()
615
616 if math.sin(wirePhi) >= 0:
617 continue
618
619 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdX), +math.cos(wirePhi))
620 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdY), +math.sin(wirePhi))
621 const.add(self.get_label(layer, 0, 0), 0.)
622 consts.append(const)
623
624 for layer in self.hemisphere:
625 const = Constraint()
626 # sum of wire rotations (BWD) in layer
627 for wire in range(0, self.wires_in_layer[layer]):
628
629 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getBackwardPos3D().phi()
630
631 if math.sin(wirePhi) < 0:
632 continue
633
634 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdX), +math.cos(wirePhi))
635 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdY), +math.sin(wirePhi))
636 const.add(self.get_label(layer, 0, 0), 0.)
637 consts.append(const)
638
639 for layer in self.hemisphere:
640 const = Constraint()
641 # sum of wire rotations (FWD) in layer
642 for wire in range(0, self.wires_in_layer[layer]):
643
644 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getForwardPos3D().phi()
645
646 if math.sin(wirePhi) >= 0:
647 continue
648
649 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdX), +math.cos(wirePhi))
650 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdY), +math.sin(wirePhi))
651 const.add(self.get_label(layer, 0, 0), 0.)
652 consts.append(const)
653
654 for layer in self.hemisphere:
655 const = Constraint()
656 # sum of wire rotations (FWD) in layer
657 for wire in range(0, self.wires_in_layer[layer]):
658
659 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getForwardPos3D().phi()
660
661 if math.sin(wirePhi) < 0:
662 continue
663
664 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdX), +math.cos(wirePhi))
665 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdY), +math.sin(wirePhi))
666 const.add(self.get_label(layer, 0, 0), 0.)
667 consts.append(const)
668
669 if self.cdc_radius:
670 const = Constraint()
671 for layer in layers:
672 # sum of wire rotations (BWD) in layer
673 for wire in range(0, self.wires_in_layer[layer]):
674 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getBackwardPos3D().phi()
675 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdX), +math.cos(wirePhi))
676 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireBwdY), +math.sin(wirePhi))
677 consts.append(const)
678
679 const = Constraint()
680 for layer in layers:
681 # sum of wire rotations (FWD) in layer
682 for wire in range(0, self.wires_in_layer[layer]):
683 wirePhi = Belle2.TrackFindingCDC.CDCWire.getInstance(Belle2.WireID(layer, wire)).getForwardPos3D().phi()
684 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdX), +math.cos(wirePhi))
685 const.add(self.get_label(layer, wire, Belle2.CDCAlignment.wireFwdY), +math.sin(wirePhi))
686 consts.append(const)
687
688 return consts
689
690
691# ------------ Main: Generate some constraint files with default config (no time-dependence, default global tags) ------
692if __name__ == '__main__':
693 ## 6D CDC constraints
694 consts6 = CDCLayerConstraints('cdc-layer-constraints-6D.txt', rigid=True, z_offset=False, r_scale=False, z_scale=False)
695 ## 7D CDC constraints
696 consts7 = CDCLayerConstraints('cdc-layer-constraints-7D.txt', rigid=True, z_offset=True, r_scale=False, z_scale=False)
697 ## 10D CD constraints
698 consts10 = CDCLayerConstraints('cdc-layer-constraints-10D.txt', rigid=True, z_offset=True, r_scale=True, z_scale=True)
699 ## CDC EventT0 constraint
700 cdcT0 = CDCTimeZerosConstraint()
701 ## CDC wire constraints
702 cdcWires = CDCWireConstraints(
703 filename='cdc-wire-constraints-proc12.txt',
704 layers=None,
705 layer_rigid=True,
706 layer_radius=[53],
707 cdc_radius=True,
708 hemisphere=[55])
709
710 # phase 2
711 # timedep = [([], [(0, 0, 1002)])]
712 # early phase 3
713 # timedep = [([], [(0, 0, 1003)])]
714
715 # final detector (phase 3)
716 ## Time dependency
717 timedep = [] # [([], [(0, 0, 0)])]
718 ## Initial event
719 init_event = (0, 0, 0)
720 ## Files
721 files = generate_constraints(
722 [
723 consts6,
724 consts7,
725 consts10,
726 cdcT0,
727 cdcWires,
728 VXDHierarchyConstraints(type=1, pxd=False),
729 Constraints("my_file.txt")],
730
731 timedep=timedep,
732 global_tags=None,
733 init_event=init_event)
static unsigned short getGlobalUniqueID()
Get global unique id.
static unsigned short getGlobalUniqueID()
Get global unique id.
static std::string findFile(const std::string &path, bool silent=false)
Search for given file or directory in local or central release directory, and return absolute path if...
Class to convert to/from global labels for Millepede II to/from detector & parameter identificators.
Definition GlobalLabel.h:41
Class to identify a wire inside the CDC.
Definition WireID.h:34
__init__(self, filename='cdc-constraints.txt', rigid=True, twist=True, z_offset=False, r_scale=False, z_scale=False)
list cdc
CDC geometry hard-code -> TODO.
__init__(self, filename='cdc-T0-constraints.txt')
list wires_in_layer
CDC layer/wire configuration.
layer_rigid
6 x 56 (6/layer) constraints.
get_label(self, layer, wire, parameter)
cdc_radius
2 Constraints: Sum(dr)=0 for all wires in CDC at each end-plate -> "average CDC radius" kept same by ...
__init__(self, filename='cdc-wire-constraints.txt', layers=None, layer_rigid=True, layer_radius=None, cdc_radius=False, hemisphere=None)
__init__(self, comment="", value=0.)
configure_collector(self, collector)
__init__(self, type=2, pxd=True, svd=True)