Belle II Software  release-05-01-25
runIPVXD_noCAF.py
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 import basf2
4 import ROOT
5 from ROOT import Belle2
6 import numpy as np
7 
8 from runIPVXD_CAF import get_calibration
9 import millepede_calibration as mpc
10 
11 cal = get_calibration(dict(), [tag for tag in basf2.conditions.globaltags])
12 collector_file = 'CollectorOutput.root'
13 collector_file = mpc.collect(cal,
14  'dimuon_skim',
15  [f for f in Belle2.Environment.Instance().getInputFilesOverride()],
16  collector_file)
17 mpc.calibrate(cal, [collector_file])
18 
19 algo = cal.algorithms[0].algorithm
20 
21 
22 # Get the payloads into handy variables
23 payloads = list(algo.getPayloads())
24 vxd = None
25 for payload in payloads:
26  if payload.name == 'VXDAlignment':
27  vxd = payload.object.IsA().DynamicCast(Belle2.VXDAlignment().IsA(), payload.object, False)
28 
29 
30 # Profile plot for all determined parameters
31 profile = ROOT.TH1F(
32  "profile",
33  "correction & errors",
34  algo.result().getNoDeterminedParameters(),
35  1,
36  algo.result().getNoDeterminedParameters())
37 
38 # Define some branch variables
39 param = np.zeros(1, dtype=int)
40 value = np.zeros(1, dtype=float)
41 correction = np.zeros(1, dtype=float)
42 error = np.zeros(1, dtype=float)
43 layer = np.zeros(1, dtype=int)
44 ladder = np.zeros(1, dtype=int)
45 sensor = np.zeros(1, dtype=int)
46 segment = np.zeros(1, dtype=int)
47 x = np.zeros(1, dtype=float)
48 y = np.zeros(1, dtype=float)
49 z = np.zeros(1, dtype=float)
50 eigenweight = np.zeros(1, dtype=float)
51 
52 # Tree with VXD data
53 vxdtree = ROOT.TTree('vxd', 'VXD data')
54 vxdtree.Branch('layer', layer, 'layer/I')
55 vxdtree.Branch('ladder', ladder, 'ladder/I')
56 vxdtree.Branch('sensor', sensor, 'sensor/I')
57 vxdtree.Branch('segment', segment, 'segment/I')
58 vxdtree.Branch('param', param, 'param/I')
59 vxdtree.Branch('value', value, 'value/D')
60 vxdtree.Branch('correction', correction, 'correction/D')
61 vxdtree.Branch('error', error, 'error/D')
62 vxdtree.Branch('x', x, 'x/D')
63 vxdtree.Branch('y', y, 'y/D')
64 vxdtree.Branch('z', z, 'z/D')
65 vxdtree.Branch('eigenweight', eigenweight, 'eigenweight/D')
66 
67 corrections = Belle2.VXDAlignment()
68 errors = Belle2.VXDAlignment()
69 eigenweights = Belle2.VXDAlignment()
70 
71 # Index of determined param
72 ibin = 0
73 for ipar in range(0, algo.result().getNoParameters()):
74  if not algo.result().isParameterDetermined(ipar):
75  continue
76 
77  ibin = ibin + 1
78 
79  label = Belle2.GlobalLabel(algo.result().getParameterLabel(ipar))
80  param[0] = label.getParameterId()
81  correction[0] = algo.result().getParameterCorrection(ipar)
82  error[0] = algo.result().getParameterError(ipar)
83 
84  if (label.getUniqueId() == Belle2.VXDAlignment.getGlobalUniqueID()):
85  sid = label.getElementId()
86  pid = label.getParameterId()
87  ew = algo.result().getEigenVectorElement(0, ipar) # + algo.result().getEigenVectorElement(1, ipar)
88 
89  if vxd:
90  value[0] = vxd.get(sid, pid)
91  else:
92  value[0] = 0.
93  corrections.set(sid, pid, correction[0])
94  errors.set(sid, pid, error[0])
95  eigenweights.set(sid, pid, ew)
96 
97  layer[0] = Belle2.VxdID(sid).getLayerNumber()
98  ladder[0] = Belle2.VxdID(sid).getLadderNumber()
99  sensor[0] = Belle2.VxdID(sid).getSensorNumber()
100  segment[0] = Belle2.VxdID(sid).getSegmentNumber()
101 
102  x[0] = 0.
103  y[0] = 0.
104  z[0] = 0.
105  eigenweight[0] = ew
106  vxdtree.Fill()
107 
108  profile.SetBinContent(ibin, value[0])
109  profile.SetBinError(ibin, error[0])
110 
111 condition = 0
112 
113 if algo.result().getNoEigenPairs():
114  maxEigenValue = algo.result().getEigenNumber(algo.result().getNoEigenPairs() - 1)
115  minEigenValue = algo.result().getEigenNumber(0)
116  condition = maxEigenValue / minEigenValue
117 
118 if condition:
119  print("Condition number of the matrix: ", condition)
120 
121 diagfile = ROOT.TFile('MillepedeJobDiagnostics.root', 'recreate')
122 diagfile.cd()
123 if vxd:
124  vxd.Write('values')
125 corrections.Write('corrections')
126 errors.Write('errors')
127 eigenweights.Write('eigenweights')
128 vxdtree.Write('vxdtree')
129 profile.Write('profile')
130 diagfile.Close()
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::VXDAlignment::getGlobalUniqueID
static unsigned short getGlobalUniqueID()
Get global unique id.
Definition: VXDAlignment.h:57
Belle2::GlobalLabel
Class to convert to/from global labels for Millepede II to/from detector & parameter identificators.
Definition: GlobalLabel.h:51
Belle2::Environment::Instance
static Environment & Instance()
Static method to get a reference to the Environment instance.
Definition: Environment.cc:31
Belle2::VXDAlignment
VXD alignment (and maybe some calibration) parameters.
Definition: VXDAlignment.h:29