]>
Commit | Line | Data |
---|---|---|
1 | //--------------------------------------------------------------------------------- | |
2 | // The AliKFParticle class | |
3 | // . | |
4 | // @author S.Gorbunov, I.Kisel | |
5 | // @version 1.0 | |
6 | // @since 13.05.07 | |
7 | // | |
8 | // Class to reconstruct and store the decayed particle parameters. | |
9 | // The method is described in CBM-SOFT note 2007-003, | |
10 | // ``Reconstruction of decayed particles based on the Kalman filter'', | |
11 | // http://www.gsi.de/documents/DOC-2007-May-14-1.pdf | |
12 | // | |
13 | // This class is ALICE interface to general mathematics in AliKFParticleBase | |
14 | // | |
15 | // -= Copyright © ALICE HLT Group =- | |
16 | //_________________________________________________________________________________ | |
17 | ||
18 | #ifndef ALIKFPARTICLE_H | |
19 | #define ALIKFPARTICLE_H | |
20 | ||
21 | #include "AliKFParticleBase.h" | |
22 | #include "TMath.h" | |
23 | ||
24 | class AliVTrack; | |
25 | class AliVVertex; | |
26 | ||
27 | class AliKFParticle :public AliKFParticleBase | |
28 | { | |
29 | ||
30 | public: | |
31 | ||
32 | //* | |
33 | //* INITIALIZATION | |
34 | //* | |
35 | ||
36 | //* Set magnetic field for all particles | |
37 | ||
38 | static void SetField( Double_t Bz ); | |
39 | ||
40 | //* Constructor (empty) | |
41 | ||
42 | AliKFParticle():AliKFParticleBase(){ ; } | |
43 | ||
44 | //* Destructor (empty) | |
45 | ||
46 | ~AliKFParticle(){ ; } | |
47 | ||
48 | //* Construction of mother particle by its 2-3-4 daughters | |
49 | ||
50 | AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2 ); | |
51 | ||
52 | AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2, | |
53 | const AliKFParticle &d3 ); | |
54 | ||
55 | AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2, | |
56 | const AliKFParticle &d3, const AliKFParticle &d4 ); | |
57 | ||
58 | //* Initialisation from "cartesian" coordinates ( X Y Z Px Py Pz ) | |
59 | //* Parameters, covariance matrix, charge and PID hypothesis should be provided | |
60 | ||
61 | void Create( const Double_t Param[], const Double_t Cov[], Int_t Charge, Int_t PID ); | |
62 | ||
63 | //* Initialisation from ALICE track, PID hypothesis shoould be provided | |
64 | ||
65 | AliKFParticle( const AliVTrack &track, Int_t PID ); | |
66 | ||
67 | //* Initialisation from VVertex | |
68 | ||
69 | AliKFParticle( const AliVVertex &vertex ); | |
70 | ||
71 | //* Initialise covariance matrix and set current parameters to 0.0 | |
72 | ||
73 | void Initialize(); | |
74 | ||
75 | //* Set decay vertex parameters for linearisation | |
76 | ||
77 | void SetVtxGuess( Double_t x, Double_t y, Double_t z ); | |
78 | ||
79 | //* | |
80 | //* ACCESSORS | |
81 | //* | |
82 | ||
83 | //* Simple accessors | |
84 | ||
85 | Double_t GetX () const ; //* x of current position | |
86 | Double_t GetY () const ; //* y of current position | |
87 | Double_t GetZ () const ; //* z of current position | |
88 | Double_t GetPx () const ; //* x-compoment of 3-momentum | |
89 | Double_t GetPy () const ; //* y-compoment of 3-momentum | |
90 | Double_t GetPz () const ; //* z-compoment of 3-momentum | |
91 | Double_t GetE () const ; //* energy | |
92 | Double_t GetS () const ; //* decay length / momentum | |
93 | Int_t GetQ () const ; //* charge | |
94 | Double_t GetChi2 () const ; //* chi^2 | |
95 | Int_t GetNDF () const ; //* Number of Degrees of Freedom | |
96 | ||
97 | Double_t GetParameter ( int i ) const ; | |
98 | Double_t GetCovariance( int i ) const ; | |
99 | Double_t GetCovariance( int i, int j ) const ; | |
100 | ||
101 | //* Accessors with calculations, value returned w/o error flag | |
102 | ||
103 | Double_t GetP () const; //* momentum | |
104 | Double_t GetPt () const; //* transverse momentum | |
105 | Double_t GetEta () const; //* pseudorapidity | |
106 | Double_t GetPhi () const; //* phi | |
107 | Double_t GetMomentum () const; //* momentum (same as GetP() ) | |
108 | Double_t GetMass () const; //* mass | |
109 | Double_t GetDecayLength () const; //* decay length | |
110 | Double_t GetLifeTime () const; //* life time | |
111 | Double_t GetR () const; //* distance to the origin | |
112 | ||
113 | //* Accessors to estimated errors | |
114 | ||
115 | Double_t GetErrX () const ; //* x of current position | |
116 | Double_t GetErrY () const ; //* y of current position | |
117 | Double_t GetErrZ () const ; //* z of current position | |
118 | Double_t GetErrPx () const ; //* x-compoment of 3-momentum | |
119 | Double_t GetErrPy () const ; //* y-compoment of 3-momentum | |
120 | Double_t GetErrPz () const ; //* z-compoment of 3-momentum | |
121 | Double_t GetErrE () const ; //* energy | |
122 | Double_t GetErrS () const ; //* decay length / momentum | |
123 | Double_t GetErrP () const ; //* momentum | |
124 | Double_t GetErrPt () const ; //* transverse momentum | |
125 | Double_t GetErrEta () const ; //* pseudorapidity | |
126 | Double_t GetErrPhi () const ; //* phi | |
127 | Double_t GetErrMomentum () const ; //* momentum | |
128 | Double_t GetErrMass () const ; //* mass | |
129 | Double_t GetErrDecayLength () const ; //* decay length | |
130 | Double_t GetErrLifeTime () const ; //* life time | |
131 | Double_t GetErrR () const ; //* distance to the origin | |
132 | ||
133 | //* Accessors with calculations( &value, &estimated sigma ) | |
134 | //* error flag returned (0 means no error during calculations) | |
135 | ||
136 | int GetP ( Double_t &P, Double_t &SigmaP ) const ; //* momentum | |
137 | int GetPt ( Double_t &Pt, Double_t &SigmaPt ) const ; //* transverse momentum | |
138 | int GetEta ( Double_t &Eta, Double_t &SigmaEta ) const ; //* pseudorapidity | |
139 | int GetPhi ( Double_t &Phi, Double_t &SigmaPhi ) const ; //* phi | |
140 | int GetMomentum ( Double_t &P, Double_t &SigmaP ) const ; //* momentum | |
141 | int GetMass ( Double_t &M, Double_t &SigmaM ) const ; //* mass | |
142 | int GetDecayLength ( Double_t &L, Double_t &SigmaL ) const ; //* decay length | |
143 | int GetLifeTime ( Double_t &T, Double_t &SigmaT ) const ; //* life time | |
144 | int GetR ( Double_t &R, Double_t &SigmaR ) const ; //* R | |
145 | ||
146 | ||
147 | //* | |
148 | //* MODIFIERS | |
149 | //* | |
150 | ||
151 | Double_t & X () ; | |
152 | Double_t & Y () ; | |
153 | Double_t & Z () ; | |
154 | Double_t & Px () ; | |
155 | Double_t & Py () ; | |
156 | Double_t & Pz () ; | |
157 | Double_t & E () ; | |
158 | Double_t & S () ; | |
159 | Int_t & Q () ; | |
160 | Double_t & Chi2 () ; | |
161 | Int_t & NDF () ; | |
162 | ||
163 | Double_t & Parameter ( int i ) ; | |
164 | Double_t & Covariance( int i ) ; | |
165 | Double_t & Covariance( int i, int j ) ; | |
166 | Double_t * Parameters () ; | |
167 | Double_t * CovarianceMatrix() ; | |
168 | ||
169 | //* | |
170 | //* CONSTRUCTION OF THE PARTICLE BY ITS DAUGHTERS AND MOTHER | |
171 | //* USING THE KALMAN FILTER METHOD | |
172 | //* | |
173 | ||
174 | ||
175 | //* Add daughter to the particle | |
176 | ||
177 | void AddDaughter( const AliKFParticle &Daughter ); | |
178 | ||
179 | //* Add daughter via += operator: ex.{ D0; D0+=Pion; D0+= Kaon; } | |
180 | ||
181 | void operator +=( const AliKFParticle &Daughter ); | |
182 | ||
183 | //* Set production vertex | |
184 | ||
185 | void SetProductionVertex( const AliKFParticle &Vtx ); | |
186 | ||
187 | //* Set mass constraint | |
188 | ||
189 | void SetMassConstraint( Double_t Mass, Double_t SigmaMass = 0 ); | |
190 | ||
191 | //* Set no decay length for resonances | |
192 | ||
193 | void SetNoDecayLength(); | |
194 | ||
195 | //* Everything in one go | |
196 | ||
197 | void Construct( const AliKFParticle *vDaughters[], int NDaughters, | |
198 | const AliKFParticle *ProdVtx=0, Double_t Mass=-1, Bool_t IsConstrained=0 ); | |
199 | ||
200 | //* | |
201 | //* TRANSPORT | |
202 | //* | |
203 | //* ( main transportation parameter is S = SignedPath/Momentum ) | |
204 | //* ( parameters of decay & production vertices are stored locally ) | |
205 | //* | |
206 | ||
207 | //* Transport the particle to its decay vertex | |
208 | ||
209 | void TransportToDecayVertex(); | |
210 | ||
211 | //* Transport the particle to its production vertex | |
212 | ||
213 | void TransportToProductionVertex(); | |
214 | ||
215 | //* Transport the particle close to xyz[] point | |
216 | ||
217 | void TransportToPoint( const Double_t xyz[] ); | |
218 | ||
219 | //* Transport the particle close to VVertex | |
220 | ||
221 | void TransportToVertex( const AliVVertex &v ); | |
222 | ||
223 | //* Transport the particle close to another particle p | |
224 | ||
225 | void TransportToParticle( const AliKFParticle &p ); | |
226 | ||
227 | //* Transport the particle on dS parameter (SignedPath/Momentum) | |
228 | ||
229 | void TransportToDS( Double_t dS ); | |
230 | ||
231 | //* Get dS to a certain space point | |
232 | ||
233 | Double_t GetDStoPoint( const Double_t xyz[] ) const ; | |
234 | ||
235 | //* Get dS to other particle p (dSp for particle p also returned) | |
236 | ||
237 | void GetDStoParticle( const AliKFParticle &p, | |
238 | Double_t &DS, Double_t &DSp ) const ; | |
239 | ||
240 | //* Get dS to other particle p in XY-plane | |
241 | ||
242 | void GetDStoParticleXY( const AliKFParticleBase &p, | |
243 | Double_t &DS, Double_t &DSp ) const ; | |
244 | ||
245 | //* | |
246 | //* OTHER UTILITIES | |
247 | //* | |
248 | ||
249 | ||
250 | //* Calculate distance from another object [cm] | |
251 | ||
252 | Double_t GetDistanceFromVertex( const Double_t vtx[] ) const ; | |
253 | Double_t GetDistanceFromVertex( const AliKFParticle &Vtx ) const ; | |
254 | Double_t GetDistanceFromVertex( const AliVVertex &Vtx ) const ; | |
255 | Double_t GetDistanceFromParticle( const AliKFParticle &p ) const ; | |
256 | ||
257 | //* Calculate sqrt(Chi2/ndf) deviation from another object | |
258 | //* ( v = [xyz]-vertex, Cv=[Cxx,Cxy,Cyy,Cxz,Cyz,Czz]-covariance matrix ) | |
259 | ||
260 | Double_t GetDeviationFromVertex( const Double_t v[], const Double_t Cv[]=0 ) const ; | |
261 | Double_t GetDeviationFromVertex( const AliKFParticle &Vtx ) const ; | |
262 | Double_t GetDeviationFromVertex( const AliVVertex &Vtx ) const ; | |
263 | Double_t GetDeviationFromParticle( const AliKFParticle &p ) const ; | |
264 | ||
265 | //* Calculate distance from another object [cm] in XY-plane | |
266 | ||
267 | Double_t GetDistanceFromVertexXY( const Double_t vtx[] ) const ; | |
268 | Double_t GetDistanceFromVertexXY( const AliKFParticle &Vtx ) const ; | |
269 | Double_t GetDistanceFromVertexXY( const AliVVertex &Vtx ) const ; | |
270 | Double_t GetDistanceFromParticleXY( const AliKFParticle &p ) const ; | |
271 | ||
272 | //* Calculate sqrt(Chi2/ndf) deviation from another object in XY plane | |
273 | //* ( v = [xyz]-vertex, Cv=[Cxx,Cxy,Cyy,Cxz,Cyz,Czz]-covariance matrix ) | |
274 | ||
275 | Double_t GetDeviationFromVertexXY( const Double_t v[], const Double_t Cv[]=0 ) const ; | |
276 | Double_t GetDeviationFromVertexXY( const AliKFParticle &Vtx ) const ; | |
277 | Double_t GetDeviationFromVertexXY( const AliVVertex &Vtx ) const ; | |
278 | Double_t GetDeviationFromParticleXY( const AliKFParticle &p ) const ; | |
279 | ||
280 | //* Calculate opennig angle between two particles | |
281 | ||
282 | Double_t GetAngle ( const AliKFParticle &p ) const ; | |
283 | Double_t GetAngleXY( const AliKFParticle &p ) const ; | |
284 | Double_t GetAngleRZ( const AliKFParticle &p ) const ; | |
285 | ||
286 | //* Subtract the particle from the vertex | |
287 | ||
288 | void SubtractFromVertex( AliKFParticle &v ) const ; | |
289 | ||
290 | protected: | |
291 | ||
292 | //* | |
293 | //* INTERNAL STUFF | |
294 | //* | |
295 | ||
296 | //* Method to access ALICE field | |
297 | ||
298 | static Double_t GetFieldAlice(); | |
299 | ||
300 | //* Other methods required by the abstract AliKFParticleBase class | |
301 | ||
302 | void GetFieldValue( const Double_t xyz[], Double_t B[] ) const ; | |
303 | void GetDStoParticle( const AliKFParticleBase &p, Double_t &DS, Double_t &DSp )const ; | |
304 | void Transport( Double_t dS, Double_t P[], Double_t C[] ) const ; | |
305 | static void GetExternalTrackParam( const AliKFParticleBase &p, Double_t &X, Double_t &Alpha, Double_t P[5] ) ; | |
306 | ||
307 | //void GetDStoParticleALICE( const AliKFParticleBase &p, Double_t &DS, Double_t &DS1 ) const; | |
308 | ||
309 | ||
310 | private: | |
311 | ||
312 | static Double_t fgBz; //* Bz compoment of the magnetic field | |
313 | ||
314 | ClassDef( AliKFParticle, 1 ); | |
315 | ||
316 | }; | |
317 | ||
318 | ||
319 | ||
320 | //--------------------------------------------------------------------- | |
321 | // | |
322 | // Inline implementation of the AliKFParticle methods | |
323 | // | |
324 | //--------------------------------------------------------------------- | |
325 | ||
326 | ||
327 | inline void AliKFParticle::SetField( Double_t Bz ) | |
328 | { | |
329 | fgBz = Bz; | |
330 | } | |
331 | ||
332 | ||
333 | inline AliKFParticle::AliKFParticle( const AliKFParticle &d1, | |
334 | const AliKFParticle &d2 ) | |
335 | { | |
336 | AliKFParticle mother; | |
337 | mother+= d1; | |
338 | mother+= d2; | |
339 | *this = mother; | |
340 | } | |
341 | ||
342 | inline AliKFParticle::AliKFParticle( const AliKFParticle &d1, | |
343 | const AliKFParticle &d2, | |
344 | const AliKFParticle &d3 ) | |
345 | { | |
346 | AliKFParticle mother; | |
347 | mother+= d1; | |
348 | mother+= d2; | |
349 | mother+= d3; | |
350 | *this = mother; | |
351 | } | |
352 | ||
353 | inline AliKFParticle::AliKFParticle( const AliKFParticle &d1, | |
354 | const AliKFParticle &d2, | |
355 | const AliKFParticle &d3, | |
356 | const AliKFParticle &d4 ) | |
357 | { | |
358 | AliKFParticle mother; | |
359 | mother+= d1; | |
360 | mother+= d2; | |
361 | mother+= d3; | |
362 | mother+= d4; | |
363 | *this = mother; | |
364 | } | |
365 | ||
366 | ||
367 | inline void AliKFParticle::Initialize() | |
368 | { | |
369 | AliKFParticleBase::Initialize(); | |
370 | } | |
371 | ||
372 | inline void AliKFParticle::SetVtxGuess( Double_t x, Double_t y, Double_t z ) | |
373 | { | |
374 | AliKFParticleBase::SetVtxGuess(x,y,z); | |
375 | } | |
376 | ||
377 | inline Double_t AliKFParticle::GetX () const | |
378 | { | |
379 | return AliKFParticleBase::GetX(); | |
380 | } | |
381 | ||
382 | inline Double_t AliKFParticle::GetY () const | |
383 | { | |
384 | return AliKFParticleBase::GetY(); | |
385 | } | |
386 | ||
387 | inline Double_t AliKFParticle::GetZ () const | |
388 | { | |
389 | return AliKFParticleBase::GetZ(); | |
390 | } | |
391 | ||
392 | inline Double_t AliKFParticle::GetPx () const | |
393 | { | |
394 | return AliKFParticleBase::GetPx(); | |
395 | } | |
396 | ||
397 | inline Double_t AliKFParticle::GetPy () const | |
398 | { | |
399 | return AliKFParticleBase::GetPy(); | |
400 | } | |
401 | ||
402 | inline Double_t AliKFParticle::GetPz () const | |
403 | { | |
404 | return AliKFParticleBase::GetPz(); | |
405 | } | |
406 | ||
407 | inline Double_t AliKFParticle::GetE () const | |
408 | { | |
409 | return AliKFParticleBase::GetE(); | |
410 | } | |
411 | ||
412 | inline Double_t AliKFParticle::GetS () const | |
413 | { | |
414 | return AliKFParticleBase::GetS(); | |
415 | } | |
416 | ||
417 | inline Int_t AliKFParticle::GetQ () const | |
418 | { | |
419 | return AliKFParticleBase::GetQ(); | |
420 | } | |
421 | ||
422 | inline Double_t AliKFParticle::GetChi2 () const | |
423 | { | |
424 | return AliKFParticleBase::GetChi2(); | |
425 | } | |
426 | ||
427 | inline Int_t AliKFParticle::GetNDF () const | |
428 | { | |
429 | return AliKFParticleBase::GetNDF(); | |
430 | } | |
431 | ||
432 | inline Double_t AliKFParticle::GetParameter ( int i ) const | |
433 | { | |
434 | return AliKFParticleBase::GetParameter(i); | |
435 | } | |
436 | ||
437 | inline Double_t AliKFParticle::GetCovariance( int i ) const | |
438 | { | |
439 | return AliKFParticleBase::GetCovariance(i); | |
440 | } | |
441 | ||
442 | inline Double_t AliKFParticle::GetCovariance( int i, int j ) const | |
443 | { | |
444 | return AliKFParticleBase::GetCovariance(i,j); | |
445 | } | |
446 | ||
447 | ||
448 | inline Double_t AliKFParticle::GetP () const | |
449 | { | |
450 | Double_t par, err; | |
451 | if( AliKFParticleBase::GetMomentum( par, err ) ) return 0; | |
452 | else return par; | |
453 | } | |
454 | ||
455 | inline Double_t AliKFParticle::GetPt () const | |
456 | { | |
457 | Double_t par, err; | |
458 | if( AliKFParticleBase::GetPt( par, err ) ) return 0; | |
459 | else return par; | |
460 | } | |
461 | ||
462 | inline Double_t AliKFParticle::GetEta () const | |
463 | { | |
464 | Double_t par, err; | |
465 | if( AliKFParticleBase::GetEta( par, err ) ) return 0; | |
466 | else return par; | |
467 | } | |
468 | ||
469 | inline Double_t AliKFParticle::GetPhi () const | |
470 | { | |
471 | Double_t par, err; | |
472 | if( AliKFParticleBase::GetPhi( par, err ) ) return 0; | |
473 | else return par; | |
474 | } | |
475 | ||
476 | inline Double_t AliKFParticle::GetMomentum () const | |
477 | { | |
478 | Double_t par, err; | |
479 | if( AliKFParticleBase::GetMomentum( par, err ) ) return 0; | |
480 | else return par; | |
481 | } | |
482 | ||
483 | inline Double_t AliKFParticle::GetMass () const | |
484 | { | |
485 | Double_t par, err; | |
486 | if( AliKFParticleBase::GetMass( par, err ) ) return 0; | |
487 | else return par; | |
488 | } | |
489 | ||
490 | inline Double_t AliKFParticle::GetDecayLength () const | |
491 | { | |
492 | Double_t par, err; | |
493 | if( AliKFParticleBase::GetDecayLength( par, err ) ) return 0; | |
494 | else return par; | |
495 | } | |
496 | ||
497 | inline Double_t AliKFParticle::GetLifeTime () const | |
498 | { | |
499 | Double_t par, err; | |
500 | if( AliKFParticleBase::GetLifeTime( par, err ) ) return 0; | |
501 | else return par; | |
502 | } | |
503 | ||
504 | inline Double_t AliKFParticle::GetR () const | |
505 | { | |
506 | Double_t par, err; | |
507 | if( AliKFParticleBase::GetR( par, err ) ) return 0; | |
508 | else return par; | |
509 | } | |
510 | ||
511 | inline Double_t AliKFParticle::GetErrX () const | |
512 | { | |
513 | return TMath::Sqrt(TMath::Abs( GetCovariance(0,0) )); | |
514 | } | |
515 | ||
516 | inline Double_t AliKFParticle::GetErrY () const | |
517 | { | |
518 | return TMath::Sqrt(TMath::Abs( GetCovariance(1,1) )); | |
519 | } | |
520 | ||
521 | inline Double_t AliKFParticle::GetErrZ () const | |
522 | { | |
523 | return TMath::Sqrt(TMath::Abs( GetCovariance(2,2) )); | |
524 | } | |
525 | ||
526 | inline Double_t AliKFParticle::GetErrPx () const | |
527 | { | |
528 | return TMath::Sqrt(TMath::Abs( GetCovariance(3,3) )); | |
529 | } | |
530 | ||
531 | inline Double_t AliKFParticle::GetErrPy () const | |
532 | { | |
533 | return TMath::Sqrt(TMath::Abs( GetCovariance(4,4) )); | |
534 | } | |
535 | ||
536 | inline Double_t AliKFParticle::GetErrPz () const | |
537 | { | |
538 | return TMath::Sqrt(TMath::Abs( GetCovariance(5,5) )); | |
539 | } | |
540 | ||
541 | inline Double_t AliKFParticle::GetErrE () const | |
542 | { | |
543 | return TMath::Sqrt(TMath::Abs( GetCovariance(6,6) )); | |
544 | } | |
545 | ||
546 | inline Double_t AliKFParticle::GetErrS () const | |
547 | { | |
548 | return TMath::Sqrt(TMath::Abs( GetCovariance(7,7) )); | |
549 | } | |
550 | ||
551 | inline Double_t AliKFParticle::GetErrP () const | |
552 | { | |
553 | Double_t par, err; | |
554 | if( AliKFParticleBase::GetMomentum( par, err ) ) return 1.e10; | |
555 | else return err; | |
556 | } | |
557 | ||
558 | inline Double_t AliKFParticle::GetErrPt () const | |
559 | { | |
560 | Double_t par, err; | |
561 | if( AliKFParticleBase::GetPt( par, err ) ) return 1.e10; | |
562 | else return err; | |
563 | } | |
564 | ||
565 | inline Double_t AliKFParticle::GetErrEta () const | |
566 | { | |
567 | Double_t par, err; | |
568 | if( AliKFParticleBase::GetEta( par, err ) ) return 1.e10; | |
569 | else return err; | |
570 | } | |
571 | ||
572 | inline Double_t AliKFParticle::GetErrPhi () const | |
573 | { | |
574 | Double_t par, err; | |
575 | if( AliKFParticleBase::GetPhi( par, err ) ) return 1.e10; | |
576 | else return err; | |
577 | } | |
578 | ||
579 | inline Double_t AliKFParticle::GetErrMomentum () const | |
580 | { | |
581 | Double_t par, err; | |
582 | if( AliKFParticleBase::GetMomentum( par, err ) ) return 1.e10; | |
583 | else return err; | |
584 | } | |
585 | ||
586 | inline Double_t AliKFParticle::GetErrMass () const | |
587 | { | |
588 | Double_t par, err; | |
589 | if( AliKFParticleBase::GetMass( par, err ) ) return 1.e10; | |
590 | else return err; | |
591 | } | |
592 | ||
593 | inline Double_t AliKFParticle::GetErrDecayLength () const | |
594 | { | |
595 | Double_t par, err; | |
596 | if( AliKFParticleBase::GetDecayLength( par, err ) ) return 1.e10; | |
597 | else return err; | |
598 | } | |
599 | ||
600 | inline Double_t AliKFParticle::GetErrLifeTime () const | |
601 | { | |
602 | Double_t par, err; | |
603 | if( AliKFParticleBase::GetLifeTime( par, err ) ) return 1.e10; | |
604 | else return err; | |
605 | } | |
606 | ||
607 | inline Double_t AliKFParticle::GetErrR () const | |
608 | { | |
609 | Double_t par, err; | |
610 | if( AliKFParticleBase::GetR( par, err ) ) return 1.e10; | |
611 | else return err; | |
612 | } | |
613 | ||
614 | ||
615 | inline int AliKFParticle::GetP( Double_t &P, Double_t &SigmaP ) const | |
616 | { | |
617 | return AliKFParticleBase::GetMomentum( P, SigmaP ); | |
618 | } | |
619 | ||
620 | inline int AliKFParticle::GetPt( Double_t &Pt, Double_t &SigmaPt ) const | |
621 | { | |
622 | return AliKFParticleBase::GetPt( Pt, SigmaPt ); | |
623 | } | |
624 | ||
625 | inline int AliKFParticle::GetEta( Double_t &Eta, Double_t &SigmaEta ) const | |
626 | { | |
627 | return AliKFParticleBase::GetEta( Eta, SigmaEta ); | |
628 | } | |
629 | ||
630 | inline int AliKFParticle::GetPhi( Double_t &Phi, Double_t &SigmaPhi ) const | |
631 | { | |
632 | return AliKFParticleBase::GetPhi( Phi, SigmaPhi ); | |
633 | } | |
634 | ||
635 | inline int AliKFParticle::GetMomentum( Double_t &P, Double_t &SigmaP ) const | |
636 | { | |
637 | return AliKFParticleBase::GetMomentum( P, SigmaP ); | |
638 | } | |
639 | ||
640 | inline int AliKFParticle::GetMass( Double_t &M, Double_t &SigmaM ) const | |
641 | { | |
642 | return AliKFParticleBase::GetMass( M, SigmaM ); | |
643 | } | |
644 | ||
645 | inline int AliKFParticle::GetDecayLength( Double_t &L, Double_t &SigmaL ) const | |
646 | { | |
647 | return AliKFParticleBase::GetDecayLength( L, SigmaL ); | |
648 | } | |
649 | ||
650 | inline int AliKFParticle::GetLifeTime( Double_t &T, Double_t &SigmaT ) const | |
651 | { | |
652 | return AliKFParticleBase::GetLifeTime( T, SigmaT ); | |
653 | } | |
654 | ||
655 | inline int AliKFParticle::GetR( Double_t &R, Double_t &SigmaR ) const | |
656 | { | |
657 | return AliKFParticleBase::GetR( R, SigmaR ); | |
658 | } | |
659 | ||
660 | inline Double_t & AliKFParticle::X() | |
661 | { | |
662 | return AliKFParticleBase::X(); | |
663 | } | |
664 | ||
665 | inline Double_t & AliKFParticle::Y() | |
666 | { | |
667 | return AliKFParticleBase::Y(); | |
668 | } | |
669 | ||
670 | inline Double_t & AliKFParticle::Z() | |
671 | { | |
672 | return AliKFParticleBase::Z(); | |
673 | } | |
674 | ||
675 | inline Double_t & AliKFParticle::Px() | |
676 | { | |
677 | return AliKFParticleBase::Px(); | |
678 | } | |
679 | ||
680 | inline Double_t & AliKFParticle::Py() | |
681 | { | |
682 | return AliKFParticleBase::Py(); | |
683 | } | |
684 | ||
685 | inline Double_t & AliKFParticle::Pz() | |
686 | { | |
687 | return AliKFParticleBase::Pz(); | |
688 | } | |
689 | ||
690 | inline Double_t & AliKFParticle::E() | |
691 | { | |
692 | return AliKFParticleBase::E(); | |
693 | } | |
694 | ||
695 | inline Double_t & AliKFParticle::S() | |
696 | { | |
697 | return AliKFParticleBase::S(); | |
698 | } | |
699 | ||
700 | inline Int_t & AliKFParticle::Q() | |
701 | { | |
702 | return AliKFParticleBase::Q(); | |
703 | } | |
704 | ||
705 | inline Double_t & AliKFParticle::Chi2() | |
706 | { | |
707 | return AliKFParticleBase::Chi2(); | |
708 | } | |
709 | ||
710 | inline Int_t & AliKFParticle::NDF() | |
711 | { | |
712 | return AliKFParticleBase::NDF(); | |
713 | } | |
714 | ||
715 | inline Double_t & AliKFParticle::Parameter ( int i ) | |
716 | { | |
717 | return AliKFParticleBase::Parameter(i); | |
718 | } | |
719 | ||
720 | inline Double_t & AliKFParticle::Covariance( int i ) | |
721 | { | |
722 | return AliKFParticleBase::Covariance(i); | |
723 | } | |
724 | ||
725 | inline Double_t & AliKFParticle::Covariance( int i, int j ) | |
726 | { | |
727 | return AliKFParticleBase::Covariance(i,j); | |
728 | } | |
729 | ||
730 | inline Double_t * AliKFParticle::Parameters () | |
731 | { | |
732 | return fP; | |
733 | } | |
734 | ||
735 | inline Double_t * AliKFParticle::CovarianceMatrix() | |
736 | { | |
737 | return fC; | |
738 | } | |
739 | ||
740 | ||
741 | inline void AliKFParticle::operator +=( const AliKFParticle &Daughter ) | |
742 | { | |
743 | AliKFParticleBase::operator +=( Daughter ); | |
744 | } | |
745 | ||
746 | ||
747 | inline void AliKFParticle::AddDaughter( const AliKFParticle &Daughter ) | |
748 | { | |
749 | AliKFParticleBase::AddDaughter( Daughter ); | |
750 | } | |
751 | ||
752 | inline void AliKFParticle::SetProductionVertex( const AliKFParticle &Vtx ) | |
753 | { | |
754 | AliKFParticleBase::SetProductionVertex( Vtx ); | |
755 | } | |
756 | ||
757 | inline void AliKFParticle::SetMassConstraint( Double_t Mass, Double_t SigmaMass ) | |
758 | { | |
759 | AliKFParticleBase::SetMassConstraint( Mass, SigmaMass ); | |
760 | } | |
761 | ||
762 | inline void AliKFParticle::SetNoDecayLength() | |
763 | { | |
764 | AliKFParticleBase::SetNoDecayLength(); | |
765 | } | |
766 | ||
767 | inline void AliKFParticle::Construct( const AliKFParticle *vDaughters[], int NDaughters, | |
768 | const AliKFParticle *ProdVtx, Double_t Mass, Bool_t IsConstrained ) | |
769 | { | |
770 | AliKFParticleBase::Construct( ( const AliKFParticleBase**)vDaughters, NDaughters, | |
771 | ( const AliKFParticleBase*)ProdVtx, Mass, IsConstrained ); | |
772 | } | |
773 | ||
774 | inline void AliKFParticle::TransportToDecayVertex() | |
775 | { | |
776 | AliKFParticleBase::TransportToDecayVertex(); | |
777 | } | |
778 | ||
779 | inline void AliKFParticle::TransportToProductionVertex() | |
780 | { | |
781 | AliKFParticleBase::TransportToProductionVertex(); | |
782 | } | |
783 | ||
784 | inline void AliKFParticle::TransportToPoint( const Double_t xyz[] ) | |
785 | { | |
786 | TransportToDS( GetDStoPoint(xyz) ); | |
787 | } | |
788 | ||
789 | inline void AliKFParticle::TransportToVertex( const AliVVertex &v ) | |
790 | { | |
791 | TransportToPoint( AliKFParticle(v).fP ); | |
792 | } | |
793 | ||
794 | inline void AliKFParticle::TransportToParticle( const AliKFParticle &p ) | |
795 | { | |
796 | Double_t dS, dSp; | |
797 | GetDStoParticle( p, dS, dSp ); | |
798 | TransportToDS( dS ); | |
799 | } | |
800 | ||
801 | inline void AliKFParticle::TransportToDS( Double_t dS ) | |
802 | { | |
803 | AliKFParticleBase::TransportToDS( dS ); | |
804 | } | |
805 | ||
806 | inline Double_t AliKFParticle::GetDStoPoint( const Double_t xyz[] ) const | |
807 | { | |
808 | return AliKFParticleBase::GetDStoPointBz( GetFieldAlice(), xyz ); | |
809 | } | |
810 | ||
811 | ||
812 | inline void AliKFParticle::GetDStoParticle( const AliKFParticle &p, | |
813 | Double_t &DS, Double_t &DSp ) const | |
814 | { | |
815 | GetDStoParticleXY( p, DS, DSp ); | |
816 | } | |
817 | ||
818 | ||
819 | inline Double_t AliKFParticle::GetDistanceFromVertex( const Double_t vtx[] ) const | |
820 | { | |
821 | return AliKFParticleBase::GetDistanceFromVertex( vtx ); | |
822 | } | |
823 | ||
824 | inline Double_t AliKFParticle::GetDeviationFromVertex( const Double_t v[], | |
825 | const Double_t Cv[] ) const | |
826 | { | |
827 | return AliKFParticleBase::GetDeviationFromVertex( v, Cv); | |
828 | } | |
829 | ||
830 | inline Double_t AliKFParticle::GetDistanceFromVertex( const AliKFParticle &Vtx ) const | |
831 | { | |
832 | return AliKFParticleBase::GetDistanceFromVertex( Vtx ); | |
833 | } | |
834 | ||
835 | inline Double_t AliKFParticle::GetDeviationFromVertex( const AliKFParticle &Vtx ) const | |
836 | { | |
837 | return AliKFParticleBase::GetDeviationFromVertex( Vtx ); | |
838 | } | |
839 | ||
840 | inline Double_t AliKFParticle::GetDistanceFromVertex( const AliVVertex &Vtx ) const | |
841 | { | |
842 | return GetDistanceFromVertex( AliKFParticle(Vtx) ); | |
843 | } | |
844 | ||
845 | inline Double_t AliKFParticle::GetDeviationFromVertex( const AliVVertex &Vtx ) const | |
846 | { | |
847 | return GetDeviationFromVertex( AliKFParticle(Vtx) ); | |
848 | } | |
849 | ||
850 | inline Double_t AliKFParticle::GetDistanceFromParticle( const AliKFParticle &p ) const | |
851 | { | |
852 | return AliKFParticleBase::GetDistanceFromParticle( p ); | |
853 | } | |
854 | ||
855 | inline Double_t AliKFParticle::GetDeviationFromParticle( const AliKFParticle &p ) const | |
856 | { | |
857 | return AliKFParticleBase::GetDeviationFromParticle( p ); | |
858 | } | |
859 | ||
860 | inline void AliKFParticle::SubtractFromVertex( AliKFParticle &v ) const | |
861 | { | |
862 | AliKFParticleBase::SubtractFromVertex( v ); | |
863 | } | |
864 | ||
865 | inline Double_t AliKFParticle::GetFieldAlice() | |
866 | { | |
867 | return fgBz; | |
868 | } | |
869 | ||
870 | inline void AliKFParticle::GetFieldValue( const Double_t * /*xyz*/, Double_t B[] ) const | |
871 | { | |
872 | B[0] = B[1] = 0; | |
873 | B[2] = GetFieldAlice(); | |
874 | } | |
875 | ||
876 | inline void AliKFParticle::GetDStoParticle( const AliKFParticleBase &p, | |
877 | Double_t &DS, Double_t &DSp )const | |
878 | { | |
879 | GetDStoParticleXY( p, DS, DSp ); | |
880 | } | |
881 | ||
882 | inline void AliKFParticle::GetDStoParticleXY( const AliKFParticleBase &p, | |
883 | Double_t &DS, Double_t &DSp ) const | |
884 | { | |
885 | AliKFParticleBase::GetDStoParticleBz( GetFieldAlice(), p, DS, DSp ) ; | |
886 | //GetDStoParticleALICE( p, DS, DSp ) ; | |
887 | } | |
888 | ||
889 | inline void AliKFParticle::Transport( Double_t dS, Double_t P[], Double_t C[] ) const | |
890 | { | |
891 | AliKFParticleBase::TransportBz( GetFieldAlice(), dS, P, C ); | |
892 | } | |
893 | ||
894 | #endif |