bugfix: not implemented functions removed from class definition (Jochen); more fighti...
[u/mrichter/AliRoot.git] / HLT / TPCLib / tracking-ca / AliHLTTPCCATrackPar.cxx
CommitLineData
326c2d4b 1// @(#) $Id$
2//*************************************************************************
3// This file is property of and copyright by the ALICE HLT Project *
4// ALICE Experiment at CERN, All rights reserved. *
5// *
6// Primary Authors: Jochen Thaeder <thaeder@kip.uni-heidelberg.de> *
7// Ivan Kisel <kisel@kip.uni-heidelberg.de> *
8// for The ALICE HLT Project. *
9// *
10// Permission to use, copy, modify and distribute this software and its *
11// documentation strictly for non-commercial purposes is hereby granted *
12// without fee, provided that the above copyright notice appears in all *
13// copies and that both the copyright notice and this permission notice *
14// appear in the supporting documentation. The authors make no claims *
15// about the suitability of this software for any purpose. It is *
16// provided "as is" without express or implied warranty. *
17//*************************************************************************
18
19#include "AliHLTTPCCATrackPar.h"
20#include "TMath.h"
21
dc4788ec 22ClassImp(AliHLTTPCCATrackPar)
326c2d4b 23
24void AliHLTTPCCATrackPar::Init()
25{
26 //* Initialization
27
28 for( Int_t i=0; i<7; i++ ) fP[i] = 0;
29 for( Int_t i=0; i<28; i++ ) fC[i] = 0;
dc4788ec 30 fC[0] = fC[2] = fC[5] = 100.;
31 fC[9] = fC[14] = fC[20] = 100.;
326c2d4b 32 fC[27] = 10.;
33 fChi2 = 0;
34 fNDF = -5;
35}
36
37
38void AliHLTTPCCATrackPar::Normalize( Double_t Direction[3] )
39{
dc4788ec 40 //* Normalize the track directions
41 //* (Px,Py,Pz,qp)-> (Px,Py,Pz,qp)/sqrt(Px²+Py²+Pz²)
42 //*
326c2d4b 43
44 Double_t p2 = fP[3]*fP[3] + fP[4]*fP[4] + fP[5]*fP[5];
45 if( p2<1.e-4 ) return;
46 Double_t a2 = 1./p2;
47 Double_t a = sqrt(a2);
48
49 if( Direction && ( fP[3]*Direction[0] + fP[4]*Direction[1] + fP[5]*Direction[2] < 0 ) ) a = -a;
50
51 Double_t ex = fP[3]*a, ey = fP[4]*a, ez = fP[5]*a, qp = fP[6]*a;
52
53 fP[3] = ex;
54 fP[4] = ey;
55 fP[5] = ez;
56 fP[6] = qp;
57
58 Double_t
59 h0 = fC[ 6]*ex + fC[10]*ey + fC[15]*ez,
60 h1 = fC[ 7]*ex + fC[11]*ey + fC[16]*ez,
61 h2 = fC[ 8]*ex + fC[12]*ey + fC[17]*ez,
62 h3 = fC[ 9]*ex + fC[13]*ey + fC[18]*ez,
63 h4 = fC[13]*ex + fC[14]*ey + fC[19]*ez,
64 h5 = fC[18]*ex + fC[19]*ey + fC[20]*ez,
65 h6 = fC[24]*ex + fC[25]*ey + fC[26]*ez,
66 d = h3*ex + h4*ey + h5*ez,
67 hh = h6 - qp*d ;
68
69 fC[ 6]= a*(fC[ 6] -ex*h0); fC[ 7]= a*(fC[ 7] -ex*h1); fC[ 8]= a*(fC[ 8] -ex*h2);
70 fC[10]= a*(fC[10] -ey*h0); fC[11]= a*(fC[11] -ey*h1); fC[12]= a*(fC[12] -ey*h2);
71 fC[15]= a*(fC[15] -ez*h0); fC[16]= a*(fC[16] -ez*h1); fC[17]= a*(fC[17] -ez*h2);
72 fC[21]= a*(fC[21] -qp*h0); fC[22]= a*(fC[22] -qp*h1); fC[23]= a*(fC[23] -qp*h2);
73
74 fC[ 9]= a2*( fC[ 9] -h3*ex -h3*ex + d*ex*ex );
75 fC[13]= a2*( fC[13] -h4*ex -h3*ey + d*ey*ex );
76 fC[14]= a2*( fC[14] -h4*ey -h4*ey + d*ey*ey );
77
78 fC[18]= a2*( fC[18] -h5*ex -h3*ez + d*ez*ex );
79 fC[19]= a2*( fC[19] -h5*ey -h4*ez + d*ez*ey );
80 fC[20]= a2*( fC[20] -h5*ez -h5*ez + d*ez*ez );
81
82 fC[24]= a2*( fC[24] -qp*h3 - hh*ex );
83 fC[25]= a2*( fC[25] -qp*h4 - hh*ey );
84 fC[26]= a2*( fC[26] -qp*h5 - hh*ez );
85 fC[27]= a2*( fC[27] -qp*h6 - hh*qp );
86}
87
88
89Double_t AliHLTTPCCATrackPar::GetDsToPointBz( Double_t Bz, const Double_t xyz[3], const Double_t *T0 ) const
90{
dc4788ec 91 //* Get dS distance to the given space point for Bz field
92 //*
93
326c2d4b 94 const Double_t kCLight = 0.000299792458;
95 Double_t bq = Bz*T0[6]*kCLight;
96 Double_t pt2 = T0[3]*T0[3] + T0[4]*T0[4];
97 if( pt2<1.e-4 ) return 0;
98 Double_t dx = xyz[0] - T0[0];
99 Double_t dy = xyz[1] - T0[1];
100 Double_t a = dx*T0[3]+dy*T0[4];
101 Double_t dS = 0;
102 if( TMath::Abs(bq)<1.e-8 ) dS = a/pt2;
103 else dS = TMath::ATan2( bq*a, pt2 + bq*(dy*T0[3] -dx*T0[4]) )/bq;
104 return dS;
105}
106
107Double_t AliHLTTPCCATrackPar::GetDsToPointBz( Double_t Bz, const Double_t xyz[3] ) const
108{
109 return GetDsToPointBz( Bz, xyz, fP );
110}
111
112void AliHLTTPCCATrackPar::TransportBz( Double_t Bz, Double_t S, Double_t *T0 )
113{
114 //* Transport the particle on dS, for Bz field
115
116 const Double_t kCLight = 0.000299792458;
117 Bz = Bz*kCLight;
118 Double_t bs= Bz*S;
119 Double_t bqs= bs*T0[6];
120 Double_t s = TMath::Sin(bqs), c = TMath::Cos(bqs);
121 Double_t sB, cB, dsB, dcB;
122 if( TMath::Abs(bqs)>1.e-10){
123 sB= s/Bz/T0[6];
124 cB= (1-c)/Bz/T0[6];
125 dsB = (c*S - sB)/T0[6];
126 dcB = (s*S-cB)/T0[6];
127 }else{
128 sB = (1. - bqs*bqs/6.)*S;
129 cB = .5*sB*bqs;
130 dsB = - T0[6]*bs*bs/3.*S;
131 dcB = .5*(sB*bs - dsB*bqs);
132 }
133
134 Double_t px = T0[3];
135 Double_t py = T0[4];
136 Double_t pz = T0[5];
137
138 Double_t d[7] = { fP[0]-T0[0], fP[1]-T0[1], fP[2]-T0[2],
139 fP[3]-T0[3], fP[4]-T0[4], fP[5]-T0[5], fP[6]-T0[6] };
140
141 T0[0] = T0[0] + sB*px + cB*py;
142 T0[1] = T0[1] - cB*px + sB*py;
143 T0[2] = T0[2] + S*pz ;
144 T0[3] = c*px + s*py;
145 T0[4] = -s*px + c*py;
146 T0[5] = T0[5];
147 T0[6] = T0[6];
148
dc4788ec 149 //* The Jacobian matrix is:
150 //* J[7][7] = { {1,0,0, sB, cB, 0, dsB*px + dcB*py },
151 //* {0,1,0, -cB, sB, 0, - dcB*px + dsB*py },
152 //* {0,0,1, 0, 0, S, 0 },
153 //* {0,0,0, c, s, 0, (-s*px + c*py)*bs },
154 //* {0,0,0, -s, c, 0, (-c*px - s*py)*bs },
155 //* {0,0,0, 0, 0, 1, 0 },
156 //* {0,0,0, 0, 0, 0, 1 } };
157 //*
158 //* below the fP=T0+J*(fP-T0) and fC=J*fC*Jt operations are performed
159 //*
160
161 Double_t h0 = dsB*px + dcB*py;
162 Double_t h1 = - dcB*px + dsB*py ;
163 Double_t h3 = (-s*px + c*py)*bs ;
164 Double_t h4 = (-c*px - s*py)*bs ;
165
166
167 fP[0] = T0[0] + d[0] + sB*d[3] + cB*d[4] + h0*d[6];
168 fP[1] = T0[1] + d[1] - cB*d[3] + sB*d[4] + h1*d[6];
169 fP[2] = T0[2] + d[2] + S*d[5];
170 fP[3] = T0[3] + c*d[3] + s*d[4] + h3*d[6];
171 fP[4] = T0[4] - s*d[3] + c*d[4] + h4*d[6];
172 fP[5] = T0[5] + d[5];
173 fP[6] = T0[6] + d[6];
326c2d4b 174
dc4788ec 175
176 Double_t
177 c00 = fC[ 0],
178 c01 = fC[ 1], c02 = fC[ 2],
179 c03 = fC[ 3], c04 = fC[ 4], c05 = fC[ 5],
180 c06 = fC[ 6], c07 = fC[ 7], c08 = fC[ 8], c09 = fC[ 9],
181 c10 = fC[10], c11 = fC[11], c12 = fC[12], c13 = fC[13], c14 = fC[14],
182 c15 = fC[15], c16 = fC[16], c17 = fC[17], c18 = fC[18], c19 = fC[19], c20 = fC[20],
183 c21 = fC[21], c22 = fC[22], c23 = fC[23], c24 = fC[24], c25 = fC[25], c26 = fC[26], c27 = fC[27];
184
185 fC[ 0] = c00 + 2*c10*cB + c14*cB*cB + 2*c21*h0 + 2*c25*cB*h0 + c27*h0*h0 + 2*c06*sB + 2*c13*cB*sB + 2*c24*h0*sB + c09*sB*sB;
186 fC[ 1] = c01 - c06*cB + c21*h1 + c10*sB + sB*(c07 - c09*cB + c24*h1 + c13*sB)
187 + cB*(c11 - c13*cB + c25*h1 + c14*sB) + h0*(c22 - c24*cB + c27*h1 + c25*sB);
188 fC[ 2] = c02 - 2*c07*cB + c09*cB*cB + 2*c22*h1 - 2*c24*cB*h1 + c27*h1*h1
189 + 2*c11*sB - 2*c13*cB*sB + 2*c25*h1*sB + c14*sB*sB;
190 fC[ 3] = c03 + c12*cB + c23*h0 + c15*S + c19*cB*S + c26*h0*S + c08*sB + c18*S*sB;
191 fC[ 4] = c04 + c16*S - cB*(c08 + c18*S) + h1*(c23 + c26*S) + (c12 + c19*S)*sB;
192 fC[ 5] = c05 + S*(2*c17 + c20*S);
193 fC[ 6] = c21*h3 + c25*cB*h3 + c27*h0*h3 + c10*s + c14*cB*s + c25*h0*s + c24*h3*sB + c13*s*sB + c*(c06 + c13*cB + c24*h0 + c09*sB);
194 fC[ 7] = c*c07 + c22*h3 + c11*s - cB*(c*c09 + c24*h3 + c13*s) +
195 h1*(c*c24 + c27*h3 + c25*s) + (c*c13 + c25*h3 + c14*s)*sB;
196 fC[ 8] = c23*h3 + c12*s + c26*h3*S + c19*s*S + c*(c08 + c18*S);
197 fC[ 9] = c*c*c09 + c27*h3*h3 + 2*c*(c24*h3 + c13*s) + s*(2*c25*h3 + c14*s);
198 fC[10] = c21*h4 + c25*cB*h4 + c27*h0*h4 - c06*s - c13*cB*s - c24*h0*s
199 + c24*h4*sB - c09*s*sB + c*(c10 + c14*cB + c25*h0 + c13*sB);
200 fC[11] = c22*h4 - c24*cB*h4 + c27*h1*h4 - c07*s + c09*cB*s - c24*h1*s +
201 c25*h4*sB - c13*s*sB + c*(c11 - c13*cB + c25*h1 + c14*sB);
202 fC[12] = c23*h4 - c08*s + c26*h4*S - c18*s*S + c*(c12 + c19*S);
203 fC[13] = c*c*c13 + c27*h3*h4 - s*(c24*h3 - c25*h4 + c13*s) + c*(c25*h3 + c24*h4 - c09*s + c14*s);
204 fC[14] = c*c*c14 + 2*c*c25*h4 + c27*h4*h4 - 2*c*c13*s - 2*c24*h4*s + c09*s*s;
205 fC[15] = c15 + c19*cB + c26*h0 + c18*sB;
206 fC[16] = c16 - c18*cB + c26*h1 + c19*sB;
207 fC[17] = c17 + c20*S;
208 fC[18] = c*c18 + c26*h3 + c19*s;
209 fC[19] = c*c19 + c26*h4 - c18*s;
210 fC[20] = c20;
211 fC[21] = c21 + c25*cB + c27*h0 + c24*sB;
212 fC[22] = c22 - c24*cB + c27*h1 + c25*sB;
213 fC[23] = c23 + c26*S;
214 fC[24] = c*c24 + c27*h3 + c25*s;
215 fC[25] = c*c25 + c27*h4 - c24*s;
216 fC[26] = c26;
217 fC[27] = c27;
326c2d4b 218}
219
220void AliHLTTPCCATrackPar::TransportBz( Double_t Bz, Double_t dS )
221{
222 //* Transport the particle on dS, for Bz field
223 TransportBz( Bz, dS, fP );
224}
225
226void AliHLTTPCCATrackPar::GetConnectionMatrix( Double_t B, const Double_t p[3], Double_t G[6], const Double_t *T0 ) const
227{
228 //* Calculate connection matrix between track and point p
229 if( !G ) return;
230 const Double_t kLight = 0.000299792458;
231 B*=kLight;
232 Double_t dx = p[0]-T0[0], dy = p[1]-T0[1], dz = p[2]-T0[2];
233 Double_t px2= T0[3]*T0[3], py2= T0[4]*T0[4], pz2= T0[5]*T0[5];
234 //Double_t B2 = B*B;
235 Double_t s2 = (dx*dx + dy*dy + dz*dz);
236 Double_t p2 = px2 + py2 + pz2;
237 if( p2>1.e-4 ) s2/=p2;
238 Double_t x = T0[3]*s2;
239 Double_t xx= px2*s2, xy= x*T0[4], xz= x*T0[5], yy= py2*s2, yz= T0[4]*T0[5]*s2;
240 //Double_t Bxy= B*xy;
241 G[ 0]= xx;
242 G[ 1]= xy; G[ 2]= yy;
243 G[ 3]= xz; G[ 4]= yz; G[ 5]= pz2*s2;
244 /*
245 C[ 0]+= xx;
246 C[ 1]+= xy; C[ 2]+= yy;
247 C[ 3]+= xz; C[ 4]+= yz; C[ 5]+= pz2*s2;
248 C[ 6]+= Bxy; C[ 7]+= B*yy; C[ 8]+= B*yz; C[ 9]+=B2*yy;
249 C[10]-= B*xx; C[11]-= Bxy; C[12]-=B*xz; C[13]-=B2*xy; C[14]+=B2*xx;
250 */
251}
252
253
2b4385d9 254void AliHLTTPCCATrackPar::Filter( const Double_t m[3], const Double_t V[6], const Double_t G[6] )
326c2d4b 255{
dc4788ec 256 //* Add the measurement m to the track using the Kalman Filter mathematics
257 //* m[3] is the measurement
258 //* V[6] is the low-triangular covariance matrix
259 //* G[6] is the track->measurement "connection matrix", additional to V[]
260 //*
261
326c2d4b 262 Double_t
263 c00 = fC[ 0],
264 c10 = fC[ 1], c11 = fC[ 2],
265 c20 = fC[ 3], c21 = fC[ 4], c22 = fC[ 5],
266 c30 = fC[ 6], c31 = fC[ 7], c32 = fC[ 8],
267 c40 = fC[10], c41 = fC[11], c42 = fC[12],
268 c50 = fC[15], c51 = fC[16], c52 = fC[17],
269 c60 = fC[21], c61 = fC[22], c62 = fC[23];
270
dc4788ec 271 Double_t
326c2d4b 272 z0 = m[0]-fP[0],
273 z1 = m[1]-fP[1],
274 z2 = m[2]-fP[2];
275
dc4788ec 276 Double_t mS[6] = { c00+V[0]+G[0], c10+V[1]+G[1], c11+V[2]+G[2],
277 c20+V[3]+G[3], c21+V[4]+G[4], c22+V[5]+G[5] };
326c2d4b 278 Double_t mSi[6];
279 mSi[0] = mS[4]*mS[4] - mS[2]*mS[5];
280 mSi[1] = mS[1]*mS[5] - mS[3]*mS[4];
281 mSi[3] = mS[2]*mS[3] - mS[1]*mS[4];
282 Double_t det = (mS[0]*mSi[0] + mS[1]*mSi[1] + mS[3]*mSi[3]);
283 if( TMath::Abs(det)<1.e-10 ) return;
284 det = 1./det;
285 mSi[0] *= det;
286 mSi[1] *= det;
287 mSi[3] *= det;
288 mSi[2] = ( mS[3]*mS[3] - mS[0]*mS[5] )*det;
289 mSi[4] = ( mS[0]*mS[4] - mS[1]*mS[3] )*det;
290 mSi[5] = ( mS[1]*mS[1] - mS[0]*mS[2] )*det;
291
292 fNDF += 2;
293 fChi2 += ( +(mSi[0]*z0 + mSi[1]*z1 + mSi[3]*z2)*z0
294 +(mSi[1]*z0 + mSi[2]*z1 + mSi[4]*z2)*z1
295 +(mSi[3]*z0 + mSi[4]*z1 + mSi[5]*z2)*z2 );
296
297 Double_t k0, k1, k2 ; // k = CHtS
298
299 k0 = c00*mSi[0] + c10*mSi[1] + c20*mSi[3];
300 k1 = c00*mSi[1] + c10*mSi[2] + c20*mSi[4];
301 k2 = c00*mSi[3] + c10*mSi[4] + c20*mSi[5];
302
303 fP[ 0]+= k0*z0 + k1*z1 + k2*z2 ;
304 fC[ 0]-= k0*c00 + k1*c10 + k2*c20;
305
306 k0 = c10*mSi[0] + c11*mSi[1] + c21*mSi[3];
307 k1 = c10*mSi[1] + c11*mSi[2] + c21*mSi[4];
308 k2 = c10*mSi[3] + c11*mSi[4] + c21*mSi[5];
309
310 fP[ 1]+= k0*z0 + k1*z1 + k2*z2 ;
311 fC[ 1]-= k0*c00 + k1*c10 + k2*c20;
312 fC[ 2]-= k0*c10 + k1*c11 + k2*c21;
313
314 k0 = c20*mSi[0] + c21*mSi[1] + c22*mSi[3];
315 k1 = c20*mSi[1] + c21*mSi[2] + c22*mSi[4];
316 k2 = c20*mSi[3] + c21*mSi[4] + c22*mSi[5];
317
318 fP[ 2]+= k0*z0 + k1*z1 + k2*z2 ;
319 fC[ 3]-= k0*c00 + k1*c10 + k2*c20;
320 fC[ 4]-= k0*c10 + k1*c11 + k2*c21;
321 fC[ 5]-= k0*c20 + k1*c21 + k2*c22;
322
323 k0 = c30*mSi[0] + c31*mSi[1] + c32*mSi[3];
324 k1 = c30*mSi[1] + c31*mSi[2] + c32*mSi[4];
325 k2 = c30*mSi[3] + c31*mSi[4] + c32*mSi[5];
326
327 fP[ 3]+= k0*z0 + k1*z1 + k2*z2 ;
328 fC[ 6]-= k0*c00 + k1*c10 + k2*c20;
329 fC[ 7]-= k0*c10 + k1*c11 + k2*c21;
330 fC[ 8]-= k0*c20 + k1*c21 + k2*c22;
331 fC[ 9]-= k0*c30 + k1*c31 + k2*c32;
332
333 k0 = c40*mSi[0] + c41*mSi[1] + c42*mSi[3];
334 k1 = c40*mSi[1] + c41*mSi[2] + c42*mSi[4];
335 k2 = c40*mSi[3] + c41*mSi[4] + c42*mSi[5];
336
337 fP[ 4]+= k0*z0 + k1*z1 + k2*z2 ;
338 fC[10]-= k0*c00 + k1*c10 + k2*c20;
339 fC[11]-= k0*c10 + k1*c11 + k2*c21;
340 fC[12]-= k0*c20 + k1*c21 + k2*c22;
341 fC[13]-= k0*c30 + k1*c31 + k2*c32;
342 fC[14]-= k0*c40 + k1*c41 + k2*c42;
343
344 k0 = c50*mSi[0] + c51*mSi[1] + c52*mSi[3];
345 k1 = c50*mSi[1] + c51*mSi[2] + c52*mSi[4];
346 k2 = c50*mSi[3] + c51*mSi[4] + c52*mSi[5];
347
348 fP[ 5]+= k0*z0 + k1*z1 + k2*z2 ;
349 fC[15]-= k0*c00 + k1*c10 + k2*c20;
350 fC[16]-= k0*c10 + k1*c11 + k2*c21;
351 fC[17]-= k0*c20 + k1*c21 + k2*c22;
352 fC[18]-= k0*c30 + k1*c31 + k2*c32;
353 fC[19]-= k0*c40 + k1*c41 + k2*c42;
354 fC[20]-= k0*c50 + k1*c51 + k2*c52;
355
356 k0 = c60*mSi[0] + c61*mSi[1] + c62*mSi[3];
357 k1 = c60*mSi[1] + c61*mSi[2] + c62*mSi[4];
358 k2 = c60*mSi[3] + c61*mSi[4] + c62*mSi[5];
359
360 fP[ 6]+= k0*z0 + k1*z1 + k2*z2 ;
361 fC[21]-= k0*c00 + k1*c10 + k2*c20;
362 fC[22]-= k0*c10 + k1*c11 + k2*c21;
363 fC[23]-= k0*c20 + k1*c21 + k2*c22;
364 fC[24]-= k0*c30 + k1*c31 + k2*c32;
365 fC[25]-= k0*c40 + k1*c41 + k2*c42;
366 fC[26]-= k0*c50 + k1*c51 + k2*c52;
367 fC[27]-= k0*c60 + k1*c61 + k2*c62;
368}
369
370
371void AliHLTTPCCATrackPar::Rotate( Double_t alpha )
372{
dc4788ec 373 //* Rotate the track parameters on the alpha angle
374
326c2d4b 375 Double_t cA = TMath::Cos( alpha );
376 Double_t sA = TMath::Sin( alpha );
377 Double_t x= fP[0], y= fP[1], px= fP[3], py= fP[4];
378 fP[0] = x*cA + y*sA;
379 fP[1] =-x*sA + y*cA;
380 fP[2] = fP[2];
381 fP[3] = px*cA + py*sA;
382 fP[4] =-px*sA + py*cA;
383 fP[5] = fP[5];
384 fP[6] = fP[6];
385
386 Double_t mJ[7][7] = { { cA,sA, 0, 0, 0, 0, 0 },
387 {-sA,cA, 0, 0, 0, 0, 0 },
388 { 0, 0, 1, 0, 0, 0, 0 },
389 { 0, 0, 0, cA, sA, 0, 0 },
390 { 0, 0, 0,-sA, cA, 0, 0 },
391 { 0, 0, 0, 0, 0, 1, 0 },
392 { 0, 0, 0, 0, 0, 0, 1 } };
393
394 Double_t mA[7][7];
395 for( Int_t k=0,i=0; i<7; i++)
396 for( Int_t j=0; j<=i; j++, k++ ) mA[i][j] = mA[j][i] = fC[k];
397
398 Double_t mJC[7][7];
399 for( Int_t i=0; i<7; i++ )
400 for( Int_t j=0; j<7; j++ ){
401 mJC[i][j]=0;
402 for( Int_t k=0; k<7; k++ ) mJC[i][j]+=mJ[i][k]*mA[k][j];
403 }
404
405 for( Int_t k=0,i=0; i<7; i++)
406 for( Int_t j=0; j<=i; j++, k++ ){
407 fC[k] = 0;
408 for( Int_t l=0; l<7; l++ ) fC[k]+=mJC[i][l]*mJ[j][l];
409 }
410}
411
412
413void AliHLTTPCCATrackPar::ConvertTo5( Double_t alpha, Double_t T[], Double_t C[] )
dc4788ec 414const
415{
416 //* Convert the track parameterisation to {y,z,ty,tz,q/p,x}
417 //* The result is stored in T[], corresponding covariance matrix in C[]
418 //*
419 //* The method is used for debuging
420 //*
421
326c2d4b 422 AliHLTTPCCATrackPar t = *this;
423 t.Rotate(alpha);
424 Double_t
425 x= t.fP[0], y= t.fP[1], z = t.fP[2],
426 ex= t.fP[3], ey= t.fP[4], ez = t.fP[5], qp = t.fP[6];
427
428 Double_t p2 = ex*ex+ey*ey+ez*ez;
429 if( p2<1.e-4 ) p2 = 1;
430 Double_t n2 = 1./p2;
431 Double_t n = sqrt(n2);
432
433 T[5] = x;
434 T[0] = y;
435 T[1] = z;
436 T[2] = ey/ex;
437 T[3] = ez/ex;
438 T[4] = qp*n;
439
440 Double_t mJ[5][7] = { { -T[2], 1, 0, 0, 0, 0, 0 },
441 { -T[3], 0, 1, 0, 0, 0, 0 },
442 { 0, 0, 0, -T[2]/ex, 1./ex, 0, 0 },
443 { 0, 0, 0, -T[3]/ex, 0, 1./ex, 0 },
444 { 0, 0, 0, -T[4]*n2*ex, -T[4]*n2*ey, -T[4]*n2*ez, n }};
445
446 Double_t mA[7][7];
447 for( Int_t k=0,i=0; i<7; i++)
448 for( Int_t j=0; j<=i; j++, k++ ) mA[i][j] = mA[j][i] = t.fC[k];
449
450 Double_t mJC[5][7];
451 for( Int_t i=0; i<5; i++ )
452 for( Int_t j=0; j<7; j++ ){
453 mJC[i][j]=0;
454 for( Int_t k=0; k<7; k++ ) mJC[i][j]+=mJ[i][k]*mA[k][j];
455 }
456
457 for( Int_t k=0,i=0; i<5; i++)
458 for( Int_t j=0; j<=i; j++, k++ ){
459 C[k] = 0;
460 for( Int_t l=0; l<7; l++ ) C[k]+=mJC[i][l]*mJ[j][l];
461 }
462}