Fork me on GitHub

source: git/external/TrackCovariance/TrkUtil.cc@ 40e890c

Last change on this file since 40e890c was a617744, checked in by michele <michele.selvaggi@…>, 4 years ago

adding latest TrackCovariance libraries (F. Bedeschi)
  • Property mode set to 100644
File size: 5.0 KB
RevLine 
[a617744]1#include "TrkUtil.h"
2
3// Constructor
4TrkUtil::TrkUtil(Double_t Bz)
5{
6 fBz = Bz;
7}
8TrkUtil::TrkUtil()
9{
10 fBz = 0.0;
11}
12//
13// Destructor
14TrkUtil::~TrkUtil()
15{
16 fBz = 0.0;
17}
18//
19// Helix parameters from position and momentum
20// static
21TVectorD TrkUtil::XPtoPar(TVector3 x, TVector3 p, Double_t Q, Double_t Bz)
22{
23 //
24 TVectorD Par(5);
25 // Transverse parameters
26 Double_t a = -Q * Bz * cSpeed(); // Units are Tesla, GeV and meters
27 Double_t pt = p.Pt();
28 Double_t C = a / (2 * pt); // Half curvature
29 //cout << "ObsTrk::XPtoPar: fB = " << fB << ", a = " << a << ", pt = " << pt << ", C = " << C << endl;
30 Double_t r2 = x.Perp2();
31 Double_t cross = x(0) * p(1) - x(1) * p(0);
32 Double_t T = TMath::Sqrt(pt * pt - 2 * a * cross + a * a * r2);
33 Double_t phi0 = TMath::ATan2((p(1) - a * x(0)) / T, (p(0) + a * x(1)) / T); // Phi0
34 Double_t D; // Impact parameter D
35 if (pt < 10.0) D = (T - pt) / a;
36 else D = (-2 * cross + a * r2) / (T + pt);
37 //
38 Par(0) = D; // Store D
39 Par(1) = phi0; // Store phi0
40 Par(2) = C; // Store C
41 //Longitudinal parameters
42 Double_t B = C * TMath::Sqrt(TMath::Max(r2 - D * D, 0.0) / (1 + 2 * C * D));
43 Double_t st = TMath::ASin(B) / C;
44 Double_t ct = p(2) / pt;
45 Double_t z0 = x(2) - ct * st;
46 //
47 Par(3) = z0; // Store z0
48 Par(4) = ct; // Store cot(theta)
49 //
50 return Par;
51}
52// non-static
53TVectorD TrkUtil::XPtoPar(TVector3 x, TVector3 p, Double_t Q)
54{
55 //
56 TVectorD Par(5);
57 Double_t Bz = fBz;
58 Par = XPtoPar(x, p, Q, Bz);
59 //
60 return Par;
61}
62//
63TVector3 TrkUtil::ParToX(TVectorD Par)
64{
65 Double_t D = Par(0);
66 Double_t phi0 = Par(1);
67 Double_t z0 = Par(3);
68 //
69 TVector3 Xval;
70 Xval(0) = -D * TMath::Sin(phi0);
71 Xval(1) = D * TMath::Cos(phi0);
72 Xval(2) = z0;
73 //
74 return Xval;
75}
76//
77TVector3 TrkUtil::ParToP(TVectorD Par)
78{
79 Double_t C = Par(2);
80 Double_t phi0 = Par(1);
81 Double_t ct = Par(4);
82 //
83 TVector3 Pval;
84 Double_t pt = fBz * cSpeed() / TMath::Abs(2 * C);
85 Pval(0) = pt * TMath::Cos(phi0);
86 Pval(1) = pt * TMath::Sin(phi0);
87 Pval(2) = pt * ct;
88 //
89 return Pval;
90}
91//
92Double_t TrkUtil::ParToQ(TVectorD Par)
93{
94 return TMath::Sign(1.0, -Par(2));
95}
96
97//
98// Parameter conversion to ACTS format
99TVectorD TrkUtil::ParToACTS(TVectorD Par)
100{
101 TVectorD pACTS(6); // Return vector
102 //
103 Double_t b = -cSpeed() * fBz / 2.;
104 pACTS(0) = 1000 * Par(0); // D from m to mm
105 pACTS(1) = 1000 * Par(3); // z0 from m to mm
106 pACTS(2) = Par(1); // Phi0 is unchanged
107 pACTS(3) = TMath::ATan2(1.0, Par(4)); // Theta in [0, pi] range
108 pACTS(4) = Par(2) / (b * TMath::Sqrt(1 + Par(4) * Par(4))); // q/p in GeV
109 pACTS(5) = 0.0; // Time: currently undefined
110 //
111 return pACTS;
112}
113// Covariance conversion to ACTS format
114TMatrixDSym TrkUtil::CovToACTS(TVectorD Par, TMatrixDSym Cov)
115{
116 TMatrixDSym cACTS(6); cACTS.Zero();
117 Double_t b = -cSpeed() * fBz / 2.;
118 //
119 // Fill derivative matrix
120 TMatrixD A(5, 5); A.Zero();
121 Double_t ct = Par(4); // cot(theta)
122 Double_t C = Par(2); // half curvature
123 A(0, 0) = 1000.; // D-D conversion to mm
124 A(1, 2) = 1.0; // phi0-phi0
125 A(2, 4) = 1.0 / (TMath::Sqrt(1.0 + ct * ct) * b); // q/p-C
126 A(3, 1) = 1000.; // z0-z0 conversion to mm
127 A(4, 3) = -1.0 / (1.0 + ct * ct); // theta - cot(theta)
128 A(4, 4) = -C * ct / (b * pow(1.0 + ct * ct, 3.0 / 2.0)); // q/p-cot(theta)
129 //
130 TMatrixDSym Cv = Cov;
131 TMatrixD At(5, 5);
132 At.Transpose(A);
133 Cv.Similarity(At);
134 TMatrixDSub(cACTS, 0, 4, 0, 4) = Cv;
135 cACTS(5, 5) = 0.1; // Currently undefined: set to arbitrary value to avoid crashes
136 //
137 return cACTS;
138}
139//
140// Parameter conversion to ILC format
141TVectorD TrkUtil::ParToILC(TVectorD Par)
142{
143 TVectorD pILC(5); // Return vector
144 //
145 pILC(0) = Par(0) * 1.0e3; // d0 in mm
146 pILC(1) = Par(1); // phi0 is unchanged
147 pILC(2) = -2 * Par(2) * 1.0e-3; // w in mm^-1
148 pILC(3) = Par(3) * 1.0e3; // z0 in mm
149 pILC(4) = Par(4); // tan(lambda) = cot(theta)
150 //
151 return pILC;
152}
153// Covariance conversion to ILC format
154TMatrixDSym TrkUtil::CovToILC(TMatrixDSym Cov)
155{
156 TMatrixDSym cILC(5); cILC.Zero();
157 //
158 // Fill derivative matrix
159 TMatrixD A(5, 5); A.Zero();
160 //
161 A(0, 0) = 1.0e3; // D-d0 in mm
162 A(1, 1) = 1.0; // phi0-phi0
163 A(2, 2) = -2.0e-3; // w-C
164 A(3, 3) = 1.0e3; // z0-z0 conversion to mm
165 A(4, 4) = 1.0; // tan(lambda) - cot(theta)
166 //
167 TMatrixDSym Cv = Cov;
168 TMatrixD At(5, 5);
169 At.Transpose(A);
170 Cv.Similarity(At);
171 cILC = Cv;
172 //
173 return cILC;
174}
175//
176// Conversion from meters to mm
177TVectorD TrkUtil::ParToMm(TVectorD Par) // Parameter conversion
178{
179 TVectorD Pmm(5); // Return vector
180 //
181 Pmm(0) = Par(0) * 1.0e3; // d0 in mm
182 Pmm(1) = Par(1); // phi0 is unchanged
183 Pmm(2) = Par(2) * 1.0e-3; // C in mm^-1
184 Pmm(3) = Par(3) * 1.0e3; // z0 in mm
185 Pmm(4) = Par(4); // tan(lambda) = cot(theta) unchanged
186 //
187 return Pmm;
188}
189TMatrixDSym TrkUtil::CovToMm(TMatrixDSym Cov) // Covariance conversion
190{
191 TMatrixDSym Cmm(5); Cmm.Zero();
192 //
193 // Fill derivative matrix
194 TMatrixD A(5, 5); A.Zero();
195 //
196 A(0, 0) = 1.0e3; // D-d0 in mm
197 A(1, 1) = 1.0; // phi0-phi0
198 A(2, 2) = 1.0e-3; // C-C
199 A(3, 3) = 1.0e3; // z0-z0 conversion to mm
200 A(4, 4) = 1.0; // lambda - cot(theta)
201 //
202 TMatrixDSym Cv = Cov;
203 TMatrixD At(5, 5);
204 At.Transpose(A);
205 Cv.Similarity(At);
206 Cmm = Cv;
207 //
208 return Cmm;
209}
Note: See TracBrowser for help on using the repository browser.