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