New TRF by Christian with x-talk
[u/mrichter/AliRoot.git] / STEER / AliExternalTrackParam.cxx
CommitLineData
51ad6848 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id$ */
17
18///////////////////////////////////////////////////////////////////////////////
19// //
20// track parameters in "external" format //
21// //
22// The track parameters are: //
23// - local y coordinate //
24// - local z coordinate //
25// - sin of azimutal angle //
26// - tan of dip angle //
27// - charge/pt //
28// The parametrisation is given at the local x coordinate fX and the //
29// azimuthal angle fAlpha. //
30// //
31// The external parametrisation can be used to exchange track parameters //
32// between different detectors. //
33// //
34///////////////////////////////////////////////////////////////////////////////
35
a0209ae4 36#include <TMath.h>
37#include <TVector3.h>
51ad6848 38
39#include "AliExternalTrackParam.h"
40#include "AliKalmanTrack.h"
a5e407e9 41#include "AliTrackReference.h"
0243e825 42#include "AliLog.h"
51ad6848 43
44ClassImp(AliExternalTrackParam)
45
a5e407e9 46const AliMagF *AliExternalTrackParam::fgkFieldMap=0;
47Double_t AliExternalTrackParam::fgConvConst=0.;
48
51ad6848 49
50//_____________________________________________________________________________
90e48c0c 51AliExternalTrackParam::AliExternalTrackParam() :
0243e825 52 fMass(-1),
90e48c0c 53 fX(0),
a5e407e9 54 fAlpha(0),
55 fLocalConvConst(0)
51ad6848 56{
90e48c0c 57 //
58 // default constructor
59 //
51ad6848 60 for (Int_t i = 0; i < 5; i++) fParam[i] = 0;
61 for (Int_t i = 0; i < 15; i++) fCovar[i] = 0;
62}
63
64//_____________________________________________________________________________
65AliExternalTrackParam::AliExternalTrackParam(Double_t x, Double_t alpha,
66 const Double_t param[5],
90e48c0c 67 const Double_t covar[15]) :
0243e825 68 fMass(-1),
90e48c0c 69 fX(x),
a5e407e9 70 fAlpha(alpha),
71 fLocalConvConst(0)
51ad6848 72{
90e48c0c 73 //
74 // create external track parameters from given arguments
75 //
51ad6848 76 for (Int_t i = 0; i < 5; i++) fParam[i] = param[i];
77 for (Int_t i = 0; i < 15; i++) fCovar[i] = covar[i];
78}
79
90e48c0c 80//_____________________________________________________________________________
81AliExternalTrackParam::AliExternalTrackParam(const AliKalmanTrack& track) :
0243e825 82 fMass(track.GetMass()),
90e48c0c 83 fX(0),
0243e825 84 fAlpha(track.GetAlpha()),
85 fLocalConvConst(0)
51ad6848 86{
87 //
88 //
89 track.GetExternalParameters(fX,fParam);
90 track.GetExternalCovariance(fCovar);
80030ed0 91 SaveLocalConvConst();
51ad6848 92}
93
a5e407e9 94
a5e407e9 95
51ad6848 96//_____________________________________________________________________________
97const Double_t* AliExternalTrackParam::GetParameter() const
98{
99// get a pointer to the array of track parameters
100
101 return fParam;
102}
103
104//_____________________________________________________________________________
105const Double_t* AliExternalTrackParam::GetCovariance() const
106{
107// get a pointer to the array of the track parameter covariance matrix
108
109 return fCovar;
110}
111
112//_____________________________________________________________________________
113AliExternalTrackParam* AliExternalTrackParam::CreateExternalParam() const
114{
115// copy this instance
116
117 return new AliExternalTrackParam(fX, fAlpha, fParam, fCovar);
118}
119
120//_____________________________________________________________________________
121void AliExternalTrackParam::ResetCovariance(Double_t factor,
122 Bool_t clearOffDiagonal)
123{
124// reset the covariance matrix ("forget" track history)
125
126 Int_t k = 0;
127 for (Int_t i = 0; i < 5; i++) {
128 for (Int_t j = 0; j < i; j++) { // off diagonal elements
129 if (clearOffDiagonal) {
130 fCovar[k++] = 0;
131 } else {
132 fCovar[k++] *= factor;
133 }
134 }
135 fCovar[k++] *= factor; // diagonal elements
136 }
137}
138
139
140//_____________________________________________________________________________
a5e407e9 141Bool_t AliExternalTrackParam::PropagateTo(Double_t xk, Double_t x0, Double_t rho)
51ad6848 142{
a5e407e9 143 //
144 // Propagate the track parameters to the given x coordinate assuming vacuum.
145 // If length is not NULL, the change of track length is added to it.
146 //
147
148 Double_t lcc=GetLocalConvConst();
149 Double_t cur = fParam[4]/lcc;
150 Double_t x1=fX, x2=xk, dx=x2-x1;
151 Double_t f1=fParam[2], f2=f1 + cur*dx;
152 if (TMath::Abs(f2) >= 0.98) {
153 // MI change - don't propagate highly inclined tracks
154 // covariance matrix distorted
155 return kFALSE;
156 }
51ad6848 157
a5e407e9 158 // old position [SR, GSI, 17.02.2003]
159 Double_t oldX = fX, oldY = fParam[0], oldZ = fParam[1];
160 Double_t r1=sqrt(1.- f1*f1), r2=sqrt(1.- f2*f2);
161 fParam[0] += dx*(f1+f2)/(r1+r2);
162 fParam[1] += dx*(f1+f2)/(f1*r2 + f2*r1)*fParam[3];
163 fParam[2] += dx*cur;
164 // transform error matrix to the curvature
165 fCovar[10]/=lcc;
166 fCovar[11]/=lcc;
167 fCovar[12]/=lcc;
168 fCovar[13]/=lcc;
169 fCovar[14]/=lcc*lcc;
170
171 //f = F - 1
172
173 Double_t f02= dx/(r1*r1*r1);
174 Double_t f04=0.5*dx*dx/(r1*r1*r1);
175 Double_t f12= dx*fParam[3]*f1/(r1*r1*r1);
176 Double_t f14=0.5*dx*dx*fParam[3]*f1/(r1*r1*r1);
177 Double_t f13= dx/r1;
178 Double_t f24= dx;
179
180 //b = C*ft
181 Double_t b00=f02*fCovar[3] + f04*fCovar[10], b01=f12*fCovar[3] + f14*fCovar[10] + f13*fCovar[6];
182 Double_t b02=f24*fCovar[10];
183 Double_t b10=f02*fCovar[4] + f04*fCovar[11], b11=f12*fCovar[4] + f14*fCovar[11] + f13*fCovar[7];
184 Double_t b12=f24*fCovar[11];
185 Double_t b20=f02*fCovar[5] + f04*fCovar[12], b21=f12*fCovar[5] + f14*fCovar[12] + f13*fCovar[8];
186 Double_t b22=f24*fCovar[12];
187 Double_t b40=f02*fCovar[12] + f04*fCovar[14], b41=f12*fCovar[12] + f14*fCovar[14] + f13*fCovar[13];
188 Double_t b42=f24*fCovar[14];
189 Double_t b30=f02*fCovar[8] + f04*fCovar[13], b31=f12*fCovar[8] + f14*fCovar[13] + f13*fCovar[9];
190 Double_t b32=f24*fCovar[13];
191
192 //a = f*b = f*C*ft
193 Double_t a00=f02*b20+f04*b40,a01=f02*b21+f04*b41,a02=f02*b22+f04*b42;
194 Double_t a11=f12*b21+f14*b41+f13*b31,a12=f12*b22+f14*b42+f13*b32;
195 Double_t a22=f24*b42;
196
197 //F*C*Ft = C + (b + bt + a)
198 fCovar[0] += b00 + b00 + a00;
199 fCovar[1] += b10 + b01 + a01;
200 fCovar[2] += b11 + b11 + a11;
201 fCovar[3] += b20 + b02 + a02;
202 fCovar[4] += b21 + b12 + a12;
203 fCovar[5] += b22 + b22 + a22;
204 fCovar[6] += b30;
205 fCovar[7] += b31;
206 fCovar[8] += b32;
207 fCovar[10] += b40;
208 fCovar[11] += b41;
209 fCovar[12] += b42;
210
211 fX=x2;
212
213 //Change of the magnetic field *************
214 SaveLocalConvConst();
215 // transform back error matrix from curvature to the 1/pt
216 fCovar[10]*=lcc;
217 fCovar[11]*=lcc;
218 fCovar[12]*=lcc;
219 fCovar[13]*=lcc;
220 fCovar[14]*=lcc*lcc;
221
222 Double_t dist = TMath::Sqrt((fX-oldX)*(fX-oldX)+(fParam[0]-oldY)*(fParam[0]-oldY)+
223 (fParam[1]-oldZ)*(fParam[1]-oldZ));
224 if (!CorrectForMaterial(dist,x0,rho)) return 0;
225
226 // Integrated Time [SR, GSI, 17.02.2003]
227 // if (IsStartedTimeIntegral() && fX>oldX) {
228// Double_t l2 = (fX-oldX)*(fX-oldX)+(fParam[0]-oldY)*(fParam[0]-oldY)+
229// (fParam[1]-oldZ)*(fParam[1]-oldZ);
230// AddTimeStep(TMath::Sqrt(l2));
231// }
232 //
233
234 return kTRUE;
235}
236
237Bool_t AliExternalTrackParam::PropagateToDCA(Double_t xd, Double_t yd, Double_t x0, Double_t rho){
238 //
239 // Propagate the track parameters to the nearest point of given xv yv coordinate
240 //
241 Double_t a=fAlpha;
242 Double_t cs=TMath::Cos(a),sn=TMath::Sin(a);
243
244 Double_t xv= xd*cs + yd*sn;
245 Double_t yv=-xd*sn + yd*cs; // vertex position in local frame
246 //
247 Double_t c=fParam[4]/GetLocalConvConst(), snp=fParam[2];
248 //
249 Double_t x=fX, y=fParam[1];
250 Double_t tgfv=-(c*(x-xv)-snp)/(c*(y-yv) + TMath::Sqrt(1.-snp*snp));
251 Double_t fv=TMath::ATan(tgfv);
252 cs=TMath::Cos(fv); sn=TMath::Sin(fv);
253 x = xv*cs + yv*sn;
254 yv=-xv*sn + yv*cs; xv=x;
255 RotateTo(fv+a);
256 return PropagateTo(xv,x0,rho);
51ad6848 257}
258
a5e407e9 259
51ad6848 260//_____________________________________________________________________________
a5e407e9 261Bool_t AliExternalTrackParam::RotateTo(Double_t alp)
51ad6848 262{
a5e407e9 263 // Rotate the reference axis for the parametrisation to the given angle.
264 //
265 Double_t x=fX;
266 Double_t p0=fParam[0];
267 //
268 if (alp < -TMath::Pi()) alp += 2*TMath::Pi();
269 else if (alp >= TMath::Pi()) alp -= 2*TMath::Pi();
270 Double_t ca=TMath::Cos(alp-fAlpha), sa=TMath::Sin(alp-fAlpha);
271 Double_t sf=fParam[2], cf=TMath::Sqrt(1.- fParam[2]*fParam[2]);
272 // **** rotation **********************
273
274 fAlpha = alp;
275 fX = x*ca + p0*sa;
276 fParam[0]= -x*sa + p0*ca;
277 fParam[2]= sf*ca - cf*sa;
278 Double_t rr=(ca+sf/cf*sa);
279 //
280 fCovar[0] *= (ca*ca);
281 fCovar[1] *= ca;
282 fCovar[3] *= ca*rr;
283 fCovar[6] *= ca;
284 fCovar[10] *= ca;
285 fCovar[4] *= rr;
286 fCovar[5] *= rr*rr;
287 fCovar[7] *= rr;
288 fCovar[11] *= rr;
289 return kTRUE;
51ad6848 290}
291
292//_____________________________________________________________________________
a5e407e9 293Bool_t AliExternalTrackParam::CorrectForMaterial(Double_t d, Double_t x0, Double_t rho)
51ad6848 294{
a5e407e9 295 //
296 // Take into account material effects assuming:
297 // x0 - mean rad length
298 // rho - mean density
299
300 //
301 // multiple scattering
302 //
0243e825 303 if (fMass<=0) {
304 AliError("Non-positive mass");
305 return kFALSE;
306 }
a5e407e9 307 Double_t p2=(1.+ fParam[3]*fParam[3])/(fParam[4]*fParam[4]);
308 Double_t beta2=p2/(p2 + fMass*fMass);
309 Double_t theta2=14.1*14.1/(beta2*p2*1e6)*d/x0*rho;
310 //
311 fCovar[5] += theta2*(1.- fParam[2]*fParam[2])*(1. + fParam[3]*fParam[3]);
312 fCovar[9] += theta2*(1. + fParam[3]*fParam[3])*(1. + fParam[3]*fParam[3]);
313 fCovar[13] += theta2*fParam[3]*fParam[4]*(1. + fParam[3]*fParam[3]);
314 fCovar[14] += theta2*fParam[3]*fParam[4]*fParam[3]*fParam[4];
315
316 Double_t dE=0.153e-3/beta2*(log(5940*beta2/(1-beta2+1e-10)) - beta2)*d*rho;
317 fParam[4] *=(1.- TMath::Sqrt(p2+fMass*fMass)/p2*dE);
318 //
319 Double_t sigmade = 0.02*TMath::Sqrt(TMath::Abs(dE)); // energy loss fluctuation
320 Double_t sigmac2 = sigmade*sigmade*fParam[4]*fParam[4]*(p2+fMass*fMass)/(p2*p2);
321 fCovar[14] += sigmac2;
322 //
323 //
324
51ad6848 325
a5e407e9 326
327 return kTRUE;
51ad6848 328}
329
330//_____________________________________________________________________________
a5e407e9 331Bool_t AliExternalTrackParam::GetProlongationAt(Double_t xk,
332 Double_t& y,
333 Double_t& z) const
51ad6848 334{
a5e407e9 335 //
336 // Get the local y and z coordinates at the given x value
337 //
338 Double_t lcc=GetLocalConvConst();
339 Double_t cur = fParam[4]/lcc;
340 Double_t x1=fX, x2=xk, dx=x2-x1;
341 Double_t f1=fParam[2], f2=f1 + cur*dx;
342 //
343 if (TMath::Abs(f2) >= 0.98) {
344 // MI change - don't propagate highly inclined tracks
345 // covariance matrix distorted
346 return kFALSE;
347 }
348 Double_t r1=sqrt(1.- f1*f1), r2=sqrt(1.- f2*f2);
349 y = fParam[0] + dx*(f1+f2)/(r1+r2);
350 z = fParam[1] + dx*(f1+f2)/(f1*r2 + f2*r1)*fParam[3];
351 return kTRUE;
51ad6848 352}
353
354//_____________________________________________________________________________
355Double_t AliExternalTrackParam::GetXAtVertex(Double_t /*x*/,
356 Double_t /*y*/) const
357{
358// Get the x coordinate at the given vertex (x,y)
359//
360// NOT IMPLEMENTED for this class
361
362 return 0;
363}
364
365
a5e407e9 366// //_____________________________________________________________________________
367// Double_t AliExternalTrackParam::GetPredictedChi2(const AliCluster* /*cluster*/)
368// {
369// // calculate the chi2 contribution of the given cluster
370// //
371// // NOT IMPLEMENTED for this class
51ad6848 372
a5e407e9 373// return -1;
374// }
51ad6848 375
a5e407e9 376// //_____________________________________________________________________________
377// Bool_t AliExternalTrackParam::Update(const AliCluster* /*cluster*/)
378// {
379// // update the track parameters using the position and error
380// // of the given cluster
381// //
382// // NOT IMPLEMENTED for this class
51ad6848 383
a5e407e9 384// return kFALSE;
385// }
51ad6848 386
387
388//_____________________________________________________________________________
389Double_t AliExternalTrackParam::SigmaPhi() const
390{
391// get the error of the azimuthal angle
392
393 return TMath::Sqrt(TMath::Abs(fCovar[5] / (1. - fParam[2]*fParam[2])));
394}
395
396//_____________________________________________________________________________
397Double_t AliExternalTrackParam::SigmaTheta() const
398{
399// get the error of the polar angle
400
401 return TMath::Sqrt(TMath::Abs(fCovar[9])) / (1. + fParam[3]*fParam[3]);
402}
403
404//_____________________________________________________________________________
405Double_t AliExternalTrackParam::SigmaPt() const
406{
407// get the error of the transversal component of the momentum
408
409 return TMath::Sqrt(fCovar[14]) / TMath::Abs(fParam[4]);
410}
411
412//_____________________________________________________________________________
413TVector3 AliExternalTrackParam::Momentum() const
414{
415// get the momentum vector
416
417 Double_t phi = TMath::ASin(fParam[2]) + fAlpha;
418 Double_t pt = 1. / TMath::Abs(fParam[4]);
419 return TVector3(pt * TMath::Cos(phi),
420 pt * TMath::Sin(phi),
421 pt * fParam[3]);
422}
423
424//_____________________________________________________________________________
425TVector3 AliExternalTrackParam::Position() const
426{
427// get the current spatial position in global coordinates
428
429 return TVector3(fX * TMath::Cos(fAlpha) - fParam[0] * TMath::Sin(fAlpha),
430 fX * TMath::Sin(fAlpha) + fParam[0] * TMath::Cos(fAlpha),
431 fParam[1]);
432}
433
434
435//_____________________________________________________________________________
436void AliExternalTrackParam::Print(Option_t* /*option*/) const
437{
438// print the parameters and the covariance matrix
439
440 printf("AliExternalTrackParam: x = %-12g alpha = %-12g\n", fX, fAlpha);
441 printf(" parameters: %12g %12g %12g %12g %12g\n",
442 fParam[0], fParam[1], fParam[2], fParam[3], fParam[4]);
443 printf(" covariance: %12g\n", fCovar[0]);
444 printf(" %12g %12g\n", fCovar[1], fCovar[2]);
445 printf(" %12g %12g %12g\n", fCovar[3], fCovar[4], fCovar[5]);
446 printf(" %12g %12g %12g %12g\n",
447 fCovar[6], fCovar[7], fCovar[8], fCovar[9]);
448 printf(" %12g %12g %12g %12g %12g\n",
449 fCovar[10], fCovar[11], fCovar[12], fCovar[13], fCovar[14]);
450}