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