]>
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 | ||
26 | ||
27 | class AliKFParticle :public AliKFParticleBase | |
28 | { | |
29 | ||
30 | public: | |
31 | ||
32 | //* | |
33 | //* INITIALIZATION | |
34 | //* | |
35 | ||
4bbc290d | 36 | //* Set magnetic field for all particles |
37 | ||
38 | static void SetField( Double_t Bz ); | |
39 | ||
f826d409 | 40 | //* Constructor (empty) |
41 | ||
4bbc290d | 42 | AliKFParticle():AliKFParticleBase(){} |
f826d409 | 43 | |
44 | //* Destructor (empty) | |
45 | ||
46 | ~AliKFParticle(){} | |
47 | ||
616ffc76 | 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 | ||
f826d409 | 58 | //* Initialisation from ALICE track, PID hypothesis can be provided |
59 | ||
616ffc76 | 60 | AliKFParticle( const AliExternalTrackParam &track, Int_t PID ); |
f826d409 | 61 | |
62 | //* Initialisation from ESD vertex | |
63 | ||
4bbc290d | 64 | AliKFParticle( const AliESDVertex &vertex ); |
f826d409 | 65 | |
66 | //* Copy vertex part to ESD vertex | |
67 | ||
68 | void CopyToESDVertex( AliESDVertex &Vtx ) const ; | |
69 | ||
70 | //* Initialise covariance matrix and set current parameters to 0.0 | |
71 | ||
72 | void Initialize(); | |
73 | ||
74 | //* Set decay vertex parameters for linearisation | |
75 | ||
76 | void SetVtxGuess( Double_t x, Double_t y, Double_t z ); | |
77 | ||
78 | //* | |
79 | //* ACCESSORS | |
80 | //* | |
81 | ||
82 | //* Simple accessors | |
83 | ||
84 | Double_t GetX () const ; //* x of current position | |
85 | Double_t GetY () const ; //* y of current position | |
86 | Double_t GetZ () const ; //* z of current position | |
87 | Double_t GetPx () const ; //* x-compoment of 3-momentum | |
88 | Double_t GetPy () const ; //* y-compoment of 3-momentum | |
89 | Double_t GetPz () const ; //* z-compoment of 3-momentum | |
90 | Double_t GetE () const ; //* energy | |
91 | Double_t GetS () const ; //* decay length / momentum | |
92 | Int_t GetQ () const ; //* charge | |
93 | Double_t GetChi2 () const ; //* chi^2 | |
94 | Int_t GetNDF () const ; //* Number of Degrees of Freedom | |
95 | ||
96 | Double_t GetParameter ( int i ) const ; | |
97 | Double_t GetCovariance( int i ) const ; | |
98 | Double_t GetCovariance( int i, int j ) const ; | |
99 | ||
100 | //* Accessors with calculations( &value, &estimated sigma ) | |
101 | //* error flag returned (0 means no error during calculations) | |
102 | ||
103 | int GetMomentum ( Double_t &P, Double_t &SigmaP ) const ; | |
104 | int GetMass ( Double_t &M, Double_t &SigmaM ) const ; | |
105 | int GetDecayLength ( Double_t &L, Double_t &SigmaL ) const ; | |
106 | int GetLifeTime ( Double_t &T, Double_t &SigmaT ) const ; | |
107 | ||
108 | //* | |
109 | //* MODIFIERS | |
110 | //* | |
111 | ||
112 | Double_t & X () ; | |
113 | Double_t & Y () ; | |
114 | Double_t & Z () ; | |
115 | Double_t & Px () ; | |
116 | Double_t & Py () ; | |
117 | Double_t & Pz () ; | |
118 | Double_t & E () ; | |
119 | Double_t & S () ; | |
120 | Int_t & Q () ; | |
121 | Double_t & Chi2 () ; | |
122 | Int_t & NDF () ; | |
123 | ||
124 | Double_t & Parameter ( int i ) ; | |
125 | Double_t & Covariance( int i ) ; | |
126 | Double_t & Covariance( int i, int j ) ; | |
127 | ||
128 | //* | |
129 | //* CONSTRUCTION OF THE PARTICLE BY ITS DAUGHTERS AND MOTHER | |
130 | //* USING THE KALMAN FILTER METHOD | |
131 | //* | |
132 | ||
133 | ||
616ffc76 | 134 | //* Add daughter to the particle |
f826d409 | 135 | |
616ffc76 | 136 | void AddDaughter( const AliKFParticle &Daughter ); |
f826d409 | 137 | |
616ffc76 | 138 | //* Add daughter via += operator: ex.{ D0; D0+=Pion; D0+= Kaon; } |
f826d409 | 139 | |
616ffc76 | 140 | void operator +=( const AliKFParticle &Daughter ); |
f826d409 | 141 | |
142 | //* Set production vertex | |
143 | ||
144 | void SetProductionVertex( const AliKFParticle &Vtx ); | |
145 | ||
146 | //* Set hard mass constraint | |
147 | ||
148 | void SetMassConstraint( Double_t Mass ); | |
149 | ||
150 | //* Everything in one go | |
151 | ||
152 | void Construct( const AliKFParticle *vDaughters[], int NDaughters, | |
153 | const AliKFParticle *ProdVtx=0, Double_t Mass=-1 ); | |
154 | ||
155 | //* | |
156 | //* TRANSPORT | |
157 | //* | |
158 | //* ( main transportation parameter is S = SignedPath/Momentum ) | |
159 | //* ( parameters of decay & production vertices are stored locally ) | |
160 | //* | |
161 | ||
162 | //* Transport the particle to its decay vertex | |
163 | ||
164 | void TransportToDecayVertex(); | |
165 | ||
166 | //* Transport the particle to its production vertex | |
167 | ||
168 | void TransportToProductionVertex(); | |
169 | ||
170 | //* Transport the particle close to xyz[] point | |
171 | ||
172 | void TransportToPoint( const Double_t xyz[] ); | |
173 | ||
174 | //* Transport the particle close to ESD vertex | |
175 | ||
176 | void TransportToVertex( const AliESDVertex &v ); | |
177 | ||
178 | //* Transport the particle close to another particle p | |
179 | ||
180 | void TransportToParticle( const AliKFParticle &p ); | |
181 | ||
182 | //* Transport the particle on dS parameter (SignedPath/Momentum) | |
183 | ||
184 | void TransportToDS( Double_t dS ); | |
185 | ||
186 | //* Get dS to a certain space point | |
187 | ||
188 | Double_t GetDStoPoint( const Double_t xyz[] ) const ; | |
189 | ||
190 | //* Get dS to other particle p (dSp for particle p also returned) | |
191 | ||
192 | void GetDStoParticle( const AliKFParticle &p, | |
193 | Double_t &DS, Double_t &DSp ) const ; | |
194 | ||
195 | //* | |
196 | //* OTHER UTILITIES | |
197 | //* | |
198 | ||
199 | ||
200 | //* Calculate distance from another object [cm] | |
201 | ||
202 | Double_t GetDistanceFromVertex( const Double_t vtx[] ) const ; | |
203 | Double_t GetDistanceFromVertex( const AliKFParticle &Vtx ) const ; | |
204 | Double_t GetDistanceFromVertex( const AliESDVertex &Vtx ) const ; | |
205 | Double_t GetDistanceFromParticle( const AliKFParticle &p ) const ; | |
206 | ||
207 | //* Calculate sqrt(Chi2/ndf) deviation from another object | |
208 | //* ( v = [xyz]-vertex, Cv=[Cxx,Cxy,Cyy,Cxz,Cyz,Czz]-covariance matrix ) | |
209 | ||
210 | Double_t GetDeviationFromVertex( const Double_t v[], const Double_t Cv[]=0 ) const ; | |
211 | Double_t GetDeviationFromVertex( const AliKFParticle &Vtx ) const ; | |
212 | Double_t GetDeviationFromVertex( const AliESDVertex &Vtx ) const ; | |
213 | Double_t GetDeviationFromParticle( const AliKFParticle &p ) const ; | |
214 | ||
215 | //* Subtract the particle from the vertex | |
216 | ||
217 | void SubtractFromVertex( AliKFParticle &v ) const ; | |
218 | void SubtractFromVertex( AliESDVertex &v ) const ; | |
219 | ||
f826d409 | 220 | protected: |
221 | ||
222 | //* | |
223 | //* INTERNAL STUFF | |
224 | //* | |
225 | ||
226 | //* Method to access ALICE field | |
227 | ||
4bbc290d | 228 | static Double_t GetFieldAlice(); |
f826d409 | 229 | |
230 | //* Other methods required by the abstract AliKFParticleBase class | |
231 | ||
232 | void GetFieldValue( const Double_t xyz[], Double_t B[] ) const ; | |
233 | void GetDStoParticle( const AliKFParticleBase &p, Double_t &DS, Double_t &DSp )const ; | |
234 | void Transport( Double_t dS, Double_t P[], Double_t C[] ) const ; | |
616ffc76 | 235 | static void GetExternalTrackParam( const AliKFParticleBase &p, Double_t &X, Double_t &Alpha, Double_t P[5] ) ; |
236 | void GetDStoParticleALICE( const AliKFParticleBase &p, Double_t &DS, Double_t &DS1 ) const; | |
f826d409 | 237 | |
4bbc290d | 238 | private: |
239 | ||
240 | static Double_t fgBz; //* Bz compoment of the magnetic field | |
f826d409 | 241 | |
616ffc76 | 242 | ClassDef( AliKFParticle, 1 ); |
f826d409 | 243 | |
244 | }; | |
245 | ||
246 | ||
247 | ||
248 | //--------------------------------------------------------------------- | |
249 | // | |
250 | // Inline implementation of the AliKFParticle methods | |
251 | // | |
252 | //--------------------------------------------------------------------- | |
253 | ||
254 | ||
4bbc290d | 255 | inline void AliKFParticle::SetField( Double_t Bz ) |
256 | { | |
616ffc76 | 257 | fgBz = -Bz;//!!! |
4bbc290d | 258 | } |
f826d409 | 259 | |
616ffc76 | 260 | |
261 | inline AliKFParticle::AliKFParticle( const AliKFParticle &d1, | |
262 | const AliKFParticle &d2 ) | |
263 | { | |
264 | AliKFParticle mother; | |
265 | mother+= d1; | |
266 | mother+= d2; | |
267 | *this = mother; | |
268 | } | |
269 | ||
270 | inline AliKFParticle::AliKFParticle( const AliKFParticle &d1, | |
271 | const AliKFParticle &d2, | |
272 | const AliKFParticle &d3 ) | |
273 | { | |
274 | AliKFParticle mother; | |
275 | mother+= d1; | |
276 | mother+= d2; | |
277 | mother+= d3; | |
278 | *this = mother; | |
279 | } | |
280 | ||
281 | inline AliKFParticle::AliKFParticle( const AliKFParticle &d1, | |
282 | const AliKFParticle &d2, | |
283 | const AliKFParticle &d3, | |
284 | const AliKFParticle &d4 ) | |
285 | { | |
286 | AliKFParticle mother; | |
287 | mother+= d1; | |
288 | mother+= d2; | |
289 | mother+= d3; | |
290 | mother+= d4; | |
291 | *this = mother; | |
292 | } | |
293 | ||
294 | ||
f826d409 | 295 | inline void AliKFParticle::Initialize() |
296 | { | |
297 | AliKFParticleBase::Initialize(); | |
298 | } | |
299 | ||
300 | inline void AliKFParticle::SetVtxGuess( Double_t x, Double_t y, Double_t z ) | |
301 | { | |
302 | AliKFParticleBase::SetVtxGuess(x,y,z); | |
303 | } | |
304 | ||
305 | inline Double_t AliKFParticle::GetX () const | |
306 | { | |
307 | return AliKFParticleBase::GetX(); | |
308 | } | |
309 | ||
310 | inline Double_t AliKFParticle::GetY () const | |
311 | { | |
312 | return AliKFParticleBase::GetY(); | |
313 | } | |
314 | ||
315 | inline Double_t AliKFParticle::GetZ () const | |
316 | { | |
317 | return AliKFParticleBase::GetZ(); | |
318 | } | |
319 | ||
320 | inline Double_t AliKFParticle::GetPx () const | |
321 | { | |
322 | return AliKFParticleBase::GetPx(); | |
323 | } | |
324 | ||
325 | inline Double_t AliKFParticle::GetPy () const | |
326 | { | |
327 | return AliKFParticleBase::GetPy(); | |
328 | } | |
329 | ||
330 | inline Double_t AliKFParticle::GetPz () const | |
331 | { | |
332 | return AliKFParticleBase::GetPz(); | |
333 | } | |
334 | ||
335 | inline Double_t AliKFParticle::GetE () const | |
336 | { | |
337 | return AliKFParticleBase::GetE(); | |
338 | } | |
339 | ||
340 | inline Double_t AliKFParticle::GetS () const | |
341 | { | |
342 | return AliKFParticleBase::GetS(); | |
343 | } | |
344 | ||
345 | inline Int_t AliKFParticle::GetQ () const | |
346 | { | |
347 | return AliKFParticleBase::GetQ(); | |
348 | } | |
349 | ||
350 | inline Double_t AliKFParticle::GetChi2 () const | |
351 | { | |
352 | return AliKFParticleBase::GetChi2(); | |
353 | } | |
354 | ||
355 | inline Int_t AliKFParticle::GetNDF () const | |
356 | { | |
357 | return AliKFParticleBase::GetNDF(); | |
358 | } | |
359 | ||
360 | inline Double_t AliKFParticle::GetParameter ( int i ) const | |
361 | { | |
362 | return AliKFParticleBase::GetParameter(i); | |
363 | } | |
364 | ||
365 | inline Double_t AliKFParticle::GetCovariance( int i ) const | |
366 | { | |
367 | return AliKFParticleBase::GetCovariance(i); | |
368 | } | |
369 | ||
370 | inline Double_t AliKFParticle::GetCovariance( int i, int j ) const | |
371 | { | |
372 | return AliKFParticleBase::GetCovariance(i,j); | |
373 | } | |
374 | ||
375 | inline int AliKFParticle::GetMomentum( Double_t &P, Double_t &SigmaP ) const | |
376 | { | |
377 | return AliKFParticleBase::GetMomentum( P, SigmaP ); | |
378 | } | |
379 | ||
380 | inline int AliKFParticle::GetMass( Double_t &M, Double_t &SigmaM ) const | |
381 | { | |
382 | return AliKFParticleBase::GetMass( M, SigmaM ); | |
383 | } | |
384 | ||
385 | inline int AliKFParticle::GetDecayLength( Double_t &L, Double_t &SigmaL ) const | |
386 | { | |
387 | return AliKFParticleBase::GetDecayLength( L, SigmaL ); | |
388 | } | |
389 | ||
390 | inline int AliKFParticle::GetLifeTime( Double_t &T, Double_t &SigmaT ) const | |
391 | { | |
392 | return AliKFParticleBase::GetLifeTime( T, SigmaT ); | |
393 | } | |
394 | ||
395 | inline Double_t & AliKFParticle::X() | |
396 | { | |
397 | return AliKFParticleBase::X(); | |
398 | } | |
399 | ||
400 | inline Double_t & AliKFParticle::Y() | |
401 | { | |
402 | return AliKFParticleBase::Y(); | |
403 | } | |
404 | ||
405 | inline Double_t & AliKFParticle::Z() | |
406 | { | |
407 | return AliKFParticleBase::Z(); | |
408 | } | |
409 | ||
410 | inline Double_t & AliKFParticle::Px() | |
411 | { | |
412 | return AliKFParticleBase::Px(); | |
413 | } | |
414 | ||
415 | inline Double_t & AliKFParticle::Py() | |
416 | { | |
417 | return AliKFParticleBase::Py(); | |
418 | } | |
419 | ||
420 | inline Double_t & AliKFParticle::Pz() | |
421 | { | |
422 | return AliKFParticleBase::Pz(); | |
423 | } | |
424 | ||
425 | inline Double_t & AliKFParticle::E() | |
426 | { | |
427 | return AliKFParticleBase::E(); | |
428 | } | |
429 | ||
430 | inline Double_t & AliKFParticle::S() | |
431 | { | |
432 | return AliKFParticleBase::S(); | |
433 | } | |
434 | ||
435 | inline Int_t & AliKFParticle::Q() | |
436 | { | |
437 | return AliKFParticleBase::Q(); | |
438 | } | |
439 | ||
440 | inline Double_t & AliKFParticle::Chi2() | |
441 | { | |
442 | return AliKFParticleBase::Chi2(); | |
443 | } | |
444 | ||
445 | inline Int_t & AliKFParticle::NDF() | |
446 | { | |
447 | return AliKFParticleBase::NDF(); | |
448 | } | |
449 | ||
450 | inline Double_t & AliKFParticle::Parameter ( int i ) | |
451 | { | |
452 | return AliKFParticleBase::Parameter(i); | |
453 | } | |
454 | ||
455 | inline Double_t & AliKFParticle::Covariance( int i ) | |
456 | { | |
457 | return AliKFParticleBase::Covariance(i); | |
458 | } | |
459 | ||
460 | inline Double_t & AliKFParticle::Covariance( int i, int j ) | |
461 | { | |
462 | return AliKFParticleBase::Covariance(i,j); | |
463 | } | |
464 | ||
465 | ||
466 | inline void AliKFParticle::operator +=( const AliKFParticle &Daughter ) | |
467 | { | |
468 | AliKFParticleBase::operator +=( Daughter ); | |
469 | } | |
470 | ||
616ffc76 | 471 | //inline AliKFParticle AliKFParticle::operator +( const AliKFParticle &Daughter ) const |
472 | //{ | |
473 | //AliKFParticle tmp; | |
474 | //tmp+= *this; | |
475 | //tmp+= Daughter; | |
476 | //return tmp; | |
477 | //} | |
f826d409 | 478 | |
479 | inline void AliKFParticle::AddDaughter( const AliKFParticle &Daughter ) | |
480 | { | |
481 | AliKFParticleBase::AddDaughter( Daughter ); | |
482 | } | |
483 | ||
484 | inline void AliKFParticle::SetProductionVertex( const AliKFParticle &Vtx ) | |
485 | { | |
486 | AliKFParticleBase::SetProductionVertex( Vtx ); | |
487 | } | |
488 | ||
489 | inline void AliKFParticle::SetMassConstraint( Double_t Mass ) | |
490 | { | |
491 | AliKFParticleBase::SetMassConstraint( Mass ); | |
492 | } | |
493 | ||
494 | inline void AliKFParticle::Construct( const AliKFParticle *vDaughters[], int NDaughters, | |
495 | const AliKFParticle *ProdVtx, Double_t Mass ) | |
496 | { | |
497 | AliKFParticleBase::Construct( ( const AliKFParticleBase**)vDaughters, NDaughters, | |
498 | ( const AliKFParticleBase*)ProdVtx, Mass ); | |
499 | } | |
500 | ||
501 | inline void AliKFParticle::TransportToDecayVertex() | |
502 | { | |
503 | AliKFParticleBase::TransportToDecayVertex(); | |
504 | } | |
505 | ||
506 | inline void AliKFParticle::TransportToProductionVertex() | |
507 | { | |
508 | AliKFParticleBase::TransportToProductionVertex(); | |
509 | } | |
510 | ||
511 | inline void AliKFParticle::TransportToPoint( const Double_t xyz[] ) | |
512 | { | |
513 | TransportToDS( GetDStoPoint(xyz) ); | |
514 | } | |
515 | ||
516 | inline void AliKFParticle::TransportToVertex( const AliESDVertex &v ) | |
517 | { | |
4bbc290d | 518 | TransportToPoint( AliKFParticle(v).fP ); |
f826d409 | 519 | } |
520 | ||
521 | inline void AliKFParticle::TransportToParticle( const AliKFParticle &p ) | |
522 | { | |
523 | Double_t dS, dSp; | |
524 | GetDStoParticle( p, dS, dSp ); | |
525 | TransportToDS( dS ); | |
526 | } | |
527 | ||
528 | inline void AliKFParticle::TransportToDS( Double_t dS ) | |
529 | { | |
530 | AliKFParticleBase::TransportToDS( dS ); | |
531 | } | |
532 | ||
533 | inline Double_t AliKFParticle::GetDStoPoint( const Double_t xyz[] ) const | |
534 | { | |
535 | return AliKFParticleBase::GetDStoPointBz( GetFieldAlice(), xyz ); | |
536 | } | |
4bbc290d | 537 | |
f826d409 | 538 | |
539 | inline void AliKFParticle::GetDStoParticle( const AliKFParticle &p, | |
4bbc290d | 540 | Double_t &DS, Double_t &DSp ) const |
f826d409 | 541 | { |
616ffc76 | 542 | GetDStoParticleALICE( p, DS, DSp ); |
543 | //AliKFParticleBase::GetDStoParticleBz( GetFieldAlice(), p, DS,DSp); | |
f826d409 | 544 | } |
545 | ||
546 | ||
547 | inline Double_t AliKFParticle::GetDistanceFromVertex( const Double_t vtx[] ) const | |
548 | { | |
549 | return AliKFParticleBase::GetDistanceFromVertex( vtx ); | |
550 | } | |
551 | ||
552 | inline Double_t AliKFParticle::GetDeviationFromVertex( const Double_t v[], | |
553 | const Double_t Cv[] ) const | |
554 | { | |
555 | return AliKFParticleBase::GetDeviationFromVertex( v, Cv); | |
556 | } | |
557 | ||
558 | inline Double_t AliKFParticle::GetDistanceFromVertex( const AliKFParticle &Vtx ) const | |
559 | { | |
560 | return AliKFParticleBase::GetDistanceFromVertex( Vtx ); | |
561 | } | |
562 | ||
563 | inline Double_t AliKFParticle::GetDeviationFromVertex( const AliKFParticle &Vtx ) const | |
564 | { | |
565 | return AliKFParticleBase::GetDeviationFromVertex( Vtx ); | |
566 | } | |
567 | ||
568 | inline Double_t AliKFParticle::GetDistanceFromVertex( const AliESDVertex &Vtx ) const | |
569 | { | |
4bbc290d | 570 | return GetDistanceFromVertex( AliKFParticle(Vtx) ); |
f826d409 | 571 | } |
572 | ||
573 | inline Double_t AliKFParticle::GetDeviationFromVertex( const AliESDVertex &Vtx ) const | |
574 | { | |
4bbc290d | 575 | return GetDeviationFromVertex( AliKFParticle(Vtx) ); |
f826d409 | 576 | } |
577 | ||
578 | inline Double_t AliKFParticle::GetDistanceFromParticle( const AliKFParticle &p ) const | |
579 | { | |
616ffc76 | 580 | return AliKFParticleBase::GetDistanceFromParticle( p ); |
f826d409 | 581 | } |
582 | ||
583 | inline Double_t AliKFParticle::GetDeviationFromParticle( const AliKFParticle &p ) const | |
584 | { | |
616ffc76 | 585 | return AliKFParticleBase::GetDeviationFromParticle( p ); |
f826d409 | 586 | } |
587 | ||
588 | inline void AliKFParticle::SubtractFromVertex( AliKFParticle &v ) const | |
589 | { | |
590 | AliKFParticleBase::SubtractFromVertex( v.fP, v.fC, v.fChi2, v.fNDF); | |
591 | } | |
592 | ||
593 | inline void AliKFParticle::SubtractFromVertex( AliESDVertex &v ) const | |
594 | { | |
4bbc290d | 595 | AliKFParticle vTmp(v); |
f826d409 | 596 | SubtractFromVertex( vTmp ); |
597 | v = AliESDVertex( vTmp.fP, vTmp.fC, vTmp.fChi2, (vTmp.fNDF +3)/2, v.GetName() ); | |
598 | } | |
599 | ||
600 | inline void AliKFParticle::CopyToESDVertex( AliESDVertex &v ) const | |
601 | { | |
602 | AliKFParticle vTmp=*this; | |
603 | v = AliESDVertex( vTmp.fP, vTmp.fC, vTmp.fChi2, (vTmp.fNDF +3)/2 ); | |
604 | } | |
605 | ||
4bbc290d | 606 | inline Double_t AliKFParticle::GetFieldAlice() |
f826d409 | 607 | { |
4bbc290d | 608 | return fgBz; |
f826d409 | 609 | } |
610 | ||
611 | inline void AliKFParticle::GetFieldValue( const Double_t * /*xyz*/, Double_t B[] ) const | |
612 | { | |
613 | B[0] = B[1] = 0; | |
614 | B[2] = GetFieldAlice(); | |
615 | } | |
616 | ||
617 | inline void AliKFParticle::GetDStoParticle( const AliKFParticleBase &p, | |
4bbc290d | 618 | Double_t &DS, Double_t &DSp )const |
f826d409 | 619 | { |
616ffc76 | 620 | GetDStoParticleALICE( p, DS, DSp ); |
621 | //AliKFParticleBase::GetDStoParticleBz( GetFieldAlice(), p, DS,DSp); | |
f826d409 | 622 | } |
623 | ||
624 | inline void AliKFParticle::Transport( Double_t dS, Double_t P[], Double_t C[] ) const | |
625 | { | |
626 | AliKFParticleBase::TransportBz( GetFieldAlice(), dS, P, C ); | |
627 | } | |
628 | ||
629 | #endif |