01-sep-2003 NvE Explicit initialisation of TObject() introduced in default constructo...
[u/mrichter/AliRoot.git] / RALICE / AliTrack.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 // $Id$
17
18 ///////////////////////////////////////////////////////////////////////////
19 // Class AliTrack
20 // Handling of the attributes of a reconstructed particle track.
21 //
22 // Coding example :
23 // ----------------
24 //
25 // Float_t a[4]={195.,1.2,-0.04,8.5};
26 // Ali4Vector pmu;
27 // pmu.SetVector(a,"car");
28 // AliTrack t1;
29 // t1.Set4Momentum(pmu);
30 //
31 // Float_t b[3]={1.2,-0.04,8.5};
32 // Ali3Vector p;
33 // p.SetVector(b,"car");
34 // AliTrack t2;
35 // t2.Set3Momentum(p);
36 // t2.SetCharge(0);
37 // t2.SetMass(1.115);
38 //
39 // t1.Data();
40 // t2.Data();
41 //
42 // Float_t pi=acos(-1.);
43 // Float_t thcms=0.2*pi; // decay theta angle in cms
44 // Float_t phicms=pi/4.; // decay theta angle in cms
45 // Float_t m1=0.938;
46 // Float_t m2=0.140;
47 // t2.Decay(m1,m2,thcms,phicms); // Track t2 decay : Lambda -> proton + pion
48 //
49 // t2.List();
50 //
51 // Int_t ndec=t2.GetNdecay();
52 // AliTrack* d1=t2.GetDecayTrack(1); // Access to decay track number 1
53 // AliTrack* d2=t2.GetDecayTrack(2); // Access to decay track number 2
54 //
55 // AliSignal s1,s2,s3,s4;
56 //
57 // .... // Code (e.g. detector readout) to fill AliSignal data
58 //
59 // AliTrack trec; // Track which will be reconstructed from signals
60 // trec.AddSignal(s1);
61 // trec.AddSignal(s3);
62 // trec.AddSignal(s4);
63 //
64 // Ali3Vector P;
65 // Float_t Q,M;
66 //
67 // ... // Code which accesses signals from trec and reconstructs
68 //        3-momentum P, charge Q, mass M etc...
69 //
70 // trec.Set3Momentum(P);
71 // trec.SetCharge(Q);
72 // trec.SetMass(M);
73 //
74 // Float_t r1[3]={1.6,-3.8,25.7};
75 // Float_t er1[3]={0.2,0.5,1.8};
76 // Float_t r2[3]={8.6,23.8,-6.7};
77 // Float_t er2[3]={0.93,1.78,0.8};
78 // AliPosition begin,end;
79 // begin.SetPosition(r1,"car");
80 // begin.SetPositionErrors(er1,"car");
81 // end.SetPosition(r2,"car");
82 // end.SetPositionErrors(er2,"car");
83 // trec.SetBeginPoint(begin);
84 // trec.SetEndPoint(end);
85 // 
86 // Note : All quantities are in GeV, GeV/c or GeV/c**2
87 //
88 //--- Author: Nick van Eijndhoven 10-jul-1997 UU-SAP Utrecht
89 //- Modified: NvE $Date$ UU-SAP Utrecht
90 ///////////////////////////////////////////////////////////////////////////
91
92 #include "AliTrack.h"
93 #include "Riostream.h"
94  
95 ClassImp(AliTrack) // Class implementation to enable ROOT I/O
96  
97 AliTrack::AliTrack() : TObject(),Ali4Vector()
98 {
99 // Default constructor
100 // All variables initialised to 0
101  Init();
102  Reset();
103 }
104 ///////////////////////////////////////////////////////////////////////////
105 void AliTrack::Init()
106 {
107 // Initialisation of pointers etc...
108  fDecays=0;
109  fSignals=0;
110  fMasses=0;
111  fDmasses=0;
112  fPmasses=0;
113  fBegin=0;
114  fEnd=0;
115  fImpactXY=0;
116  fImpactXZ=0;
117  fImpactYZ=0;
118  fClosest=0;
119  fParent=0;
120 }
121 ///////////////////////////////////////////////////////////////////////////
122 AliTrack::~AliTrack()
123 {
124 // Destructor to delete memory allocated for decay tracks array
125  if (fDecays)
126  {
127   delete fDecays;
128   fDecays=0;
129  }
130  if (fSignals)
131  {
132   fSignals->Clear();
133   delete fSignals;
134   fSignals=0;
135  }
136  if (fMasses)
137  {
138   delete fMasses;
139   fMasses=0;
140  }
141  if (fDmasses)
142  {
143   delete fDmasses;
144   fDmasses=0;
145  }
146  if (fPmasses)
147  {
148   delete fPmasses;
149   fPmasses=0;
150  }
151  if (fBegin)
152  {
153   delete fBegin;
154   fBegin=0;
155  }
156  if (fEnd)
157  {
158   delete fEnd;
159   fEnd=0;
160  }
161  if (fImpactXY)
162  {
163   delete fImpactXY;
164   fImpactXY=0;
165  }
166  if (fImpactXZ)
167  {
168   delete fImpactXZ;
169   fImpactXZ=0;
170  }
171  if (fImpactYZ)
172  {
173   delete fImpactYZ;
174   fImpactYZ=0;
175  }
176  if (fClosest)
177  {
178   delete fClosest;
179   fClosest=0;
180  }
181 }
182 ///////////////////////////////////////////////////////////////////////////
183 AliTrack::AliTrack(AliTrack& t) : TObject(t),Ali4Vector(t)
184 {
185 // Copy constructor
186  Init();
187
188  fQ=t.fQ;
189  fNdec=t.fNdec;
190  fNsig=t.fNsig;
191  fNmasses=t.fNmasses;
192  if (fNmasses)
193  {
194   fMasses=new TArrayD(*(t.fMasses));
195   fDmasses=new TArrayD(*(t.fDmasses));
196   fPmasses=new TArrayD(*(t.fPmasses));
197  }
198  if (t.fBegin) fBegin=new AliPositionObj(*(t.fBegin));
199  if (t.fEnd) fEnd=new AliPositionObj(*(t.fEnd));
200  if (t.fImpactXY) fImpactXY=new AliPositionObj(*(t.fImpactXY));
201  if (t.fImpactXZ) fImpactXZ=new AliPositionObj(*(t.fImpactXZ));
202  if (t.fImpactYZ) fImpactYZ=new AliPositionObj(*(t.fImpactYZ));
203  if (t.fClosest) fClosest=new AliPositionObj(*(t.fClosest));
204  fUserId=t.fUserId;
205  fChi2=t.fChi2;
206  fNdf=t.fNdf;
207  fCode=t.fCode;
208  fParent=t.fParent;
209
210  if (fNdec)
211  {
212   fDecays=new TObjArray(fNdec);
213   fDecays->SetOwner();
214   for (Int_t it=1; it<=fNdec; it++)
215   {
216    AliTrack* tx=t.GetDecayTrack(it);
217    fDecays->Add(new AliTrack(*tx));
218   }
219  }
220
221  if (fNsig)
222  {
223   fSignals=new TObjArray(fNsig);
224   for (Int_t is=1; is<=fNsig; is++)
225   {
226    AliSignal* sx=t.GetSignal(is);
227    fSignals->Add(sx);
228   }
229  }
230 }
231 ///////////////////////////////////////////////////////////////////////////
232 void AliTrack::Reset()
233 {
234 // Reset all variables to 0 and delete all auto-generated decay tracks.
235  fQ=0;
236  fChi2=0;
237  fNdf=0;
238  fUserId=0;
239  fCode=0;
240  fNdec=0;
241  fNsig=0;
242  fNmasses=0;
243  Double_t a[4]={0,0,0,0};
244  SetVector(a,"sph");
245  fParent=0;
246  if (fDecays)
247  {
248   delete fDecays;
249   fDecays=0;
250  }
251  if (fSignals)
252  {
253   fSignals->Clear();
254   delete fSignals;
255   fSignals=0;
256  }
257  if (fMasses)
258  {
259   delete fMasses;
260   fMasses=0;
261  }
262  if (fDmasses)
263  {
264   delete fDmasses;
265   fDmasses=0;
266  }
267  if (fPmasses)
268  {
269   delete fPmasses;
270   fPmasses=0;
271  }
272  if (fBegin)
273  {
274   delete fBegin;
275   fBegin=0;
276  }
277  if (fEnd)
278  {
279   delete fEnd;
280   fEnd=0;
281  }
282  if (fImpactXY)
283  {
284   delete fImpactXY;
285   fImpactXY=0;
286  }
287  if (fImpactXZ)
288  {
289   delete fImpactXZ;
290   fImpactXZ=0;
291  }
292  if (fImpactYZ)
293  {
294   delete fImpactYZ;
295   fImpactYZ=0;
296  }
297  if (fClosest)
298  {
299   delete fClosest;
300   fClosest=0;
301  }
302 }
303 ///////////////////////////////////////////////////////////////////////////
304 void AliTrack::Set3Momentum(Ali3Vector& p)
305 {
306 // Set the track parameters according to the 3-momentum p
307  Set3Vector(p);
308 }
309 ///////////////////////////////////////////////////////////////////////////
310 void AliTrack::Set4Momentum(Ali4Vector& p)
311 {
312 // Set the track parameters according to the 4-momentum p
313  Double_t E=p.GetScalar();
314  Double_t dE=p.GetResultError();
315  Ali3Vector pv=p.Get3Vector();
316  SetVector(E,pv);
317  SetScalarError(dE);
318 }
319 ///////////////////////////////////////////////////////////////////////////
320 void AliTrack::SetMass(Double_t m,Double_t dm)
321 {
322 // Set the particle mass
323 // The default value for the error dm is 0.
324  Double_t inv=pow(m,2);
325  Double_t dinv=fabs(2.*m*dm);
326  SetInvariant(inv,dinv);
327 }
328 ///////////////////////////////////////////////////////////////////////////
329 void AliTrack::SetCharge(Float_t q)
330 {
331 // Set the particle charge
332  fQ=q;
333 }
334 ///////////////////////////////////////////////////////////////////////////
335 void AliTrack::Data(TString f)
336 {
337 // Provide track information within the coordinate frame f
338  Double_t m=GetMass();
339  Double_t dm=GetResultError();
340  cout << " *AliTrack::Data* Id : " << fUserId << " Code : " << fCode
341       << " Mass : " << m << " error : " << dm << " Charge : " << fQ
342       << " Momentum : " << GetMomentum() << " Nmass hyp. : " << fNmasses
343       << " Ntracks : " << fNdec << " Nsignals : " << fNsig << endl;
344  for (Int_t i=0; i<fNmasses; i++)
345  {
346   cout << " Mass hypothesis " << (i+1) << " Mass : " << fMasses->At(i)
347        << " error : " << fDmasses->At(i) << " prob. : " << fPmasses->At(i)
348        << endl;
349  }
350  Ali4Vector::Data(f); 
351
352 ///////////////////////////////////////////////////////////////////////////
353 void AliTrack::List(TString f)
354 {
355 // Provide current track and decay level 1 information within coordinate frame f
356
357  Data(f); // Information of the current track
358
359  // Decay products of this track
360  AliTrack* td; 
361  for (Int_t id=1; id<=fNdec; id++)
362  {
363   td=GetDecayTrack(id);
364   if (td)
365   {
366    cout << "  ---Level 1 sec. track no. " << id << endl;
367    td->Data(f); 
368   }
369   else
370   {
371    cout << " *AliTrack::List* Error : No decay track present." << endl; 
372   }
373  }
374
375 ///////////////////////////////////////////////////////////////////////////
376 void AliTrack::ListAll(TString f)
377 {
378 // Provide complete track and decay information within the coordinate frame f
379
380  Data(f); // Information of the current track
381  if (fBegin) { cout << " Begin-point :"; fBegin->Data(f); }
382  if (fEnd)   { cout << " End-point   :"; fEnd->Data(f); }
383  for (Int_t is=1; is<=GetNsignals(); is++)
384  {
385   ((AliSignal*)GetSignal(is))->Data(f);
386  }
387
388  AliTrack* t=this;
389  Dump(t,1,f); // Information of all decay products
390 }
391 //////////////////////////////////////////////////////////////////////////
392 void AliTrack::Dump(AliTrack* t,Int_t n,TString f)
393 {
394 // Recursively provide the info of all decay levels of this track
395  AliTrack* td; 
396  for (Int_t id=1; id<=t->GetNdecay(); id++)
397  {
398   td=t->GetDecayTrack(id);
399   if (td)
400   {
401    cout << "  ---Level " << n << " sec. track no. " << id << endl;
402    td->Data(f); 
403    for (Int_t is=1; is<=td->GetNsignals(); is++)
404    {
405     ((AliSignal*)td->GetSignal(is))->Data(f);
406    }
407
408    // Go for next decay level of this decay track recursively
409    Dump(td,n+1,f);
410   }
411   else
412   {
413    cout << " *AliTrack::Dump* Error : No decay track present." << endl; 
414   }
415  }
416
417 //////////////////////////////////////////////////////////////////////////
418 Double_t AliTrack::GetMomentum()
419 {
420 // Provide the value of the track 3-momentum.
421 // The error can be obtained by invoking GetResultError() after
422 // invokation of GetMomentum().
423  Double_t norm=fV.GetNorm();
424  fDresult=fV.GetResultError();
425  return norm;
426 }
427 ///////////////////////////////////////////////////////////////////////////
428 Ali3Vector AliTrack::Get3Momentum()
429 {
430 // Provide the track 3-momentum
431  return (Ali3Vector)Get3Vector();
432 }
433 ///////////////////////////////////////////////////////////////////////////
434 Double_t AliTrack::GetMass()
435 {
436 // Provide the particle mass.
437 // The error can be obtained by invoking GetResultError() after
438 // invokation of GetMass().
439  Double_t inv=GetInvariant();
440  Double_t dinv=GetResultError();
441  Double_t dm=0;
442  if (inv >= 0)
443  {
444  Double_t m=sqrt(inv);
445  if (m) dm=dinv/(2.*m);
446  fDresult=dm;
447  return m;
448  }
449  else
450  {
451   cout << "*AliTrack::GetMass* Unphysical situation m**2 = " << inv << endl;
452   cout << " Value 0 will be returned." << endl;
453   fDresult=dm;
454   return 0;
455  }
456 }
457 ///////////////////////////////////////////////////////////////////////////
458 Float_t AliTrack::GetCharge()
459 {
460 // Provide the particle charge
461  return fQ;
462 }
463 ///////////////////////////////////////////////////////////////////////////
464 Double_t AliTrack::GetEnergy()
465 {
466 // Provide the particle's energy.
467 // The error can be obtained by invoking GetResultError() after
468 // invokation of GetEnergy().
469  Double_t E=GetScalar();
470  if (E>0)
471  {
472   return E;
473  }
474  else
475  {
476   cout << "*AliTrack::GetEnergy* Unphysical situation E = " << E << endl;
477   cout << " Value 0 will be returned." << endl;
478   return 0;
479  }
480 }
481 ///////////////////////////////////////////////////////////////////////////
482 void AliTrack::Decay(Double_t m1,Double_t m2,Double_t thcms,Double_t phicms)
483 {
484 // Perform 2-body decay of current track
485 // m1     : mass of decay product 1
486 // m2     : mass of decay product 2
487 // thcms  : cms theta decay angle (in rad.) of m1
488 // phicms : cms phi decay angle (in rad.) of m1
489  
490  fNdec=2; // it's a 2-body decay
491
492  Double_t M=GetMass();
493  
494 // Compute the 4-momenta of the decay products in the cms
495 // Note : p2=p1=pnorm for a 2-body decay
496  Double_t e1=0;
497  if (M) e1=((M*M)+(m1*m1)-(m2*m2))/(2.*M);
498  Double_t e2=0;
499  if (M) e2=((M*M)+(m2*m2)-(m1*m1))/(2.*M);
500  Double_t pnorm=(e1*e1)-(m1*m1);
501  if (pnorm>0.)
502  {
503   pnorm=sqrt(pnorm);
504  }
505  else
506  {
507   pnorm=0;
508  }
509  
510  Double_t a[3];
511  a[0]=pnorm;
512  a[1]=thcms;
513  a[2]=phicms;
514  Ali3Vector p;
515  p.SetVector(a,"sph");
516
517  Ali4Vector pprim1;
518  pprim1.SetVector(e1,p);
519  pprim1.SetInvariant(m1*m1);
520
521  Ali4Vector pprim2;
522  p*=-1;
523  pprim2.SetVector(e2,p);
524  pprim2.SetInvariant(m2*m2);
525
526  // Determine boost parameters from the parent particle
527  Double_t E=GetEnergy();
528  p=Get3Vector();
529  Ali4Vector pmu;
530  pmu.SetVector(E,p);
531
532  AliBoost q;
533  q.Set4Momentum(pmu);
534  
535  Ali4Vector p1=q.Inverse(pprim1); // Boost decay product 1
536  Ali4Vector p2=q.Inverse(pprim2); // Boost decay product 2
537  
538  // Enter the boosted data into the decay tracks array
539  if (fDecays)
540  {
541   delete fDecays;
542   fDecays=0;
543  }
544  fDecays=new TObjArray();
545  fDecays->SetOwner();
546
547  fDecays->Add(new AliTrack);
548  ((AliTrack*)fDecays->At(0))->Set4Momentum(p1);
549  ((AliTrack*)fDecays->At(0))->SetMass(m1);
550  fDecays->Add(new AliTrack);
551  ((AliTrack*)fDecays->At(1))->Set4Momentum(p2);
552  ((AliTrack*)fDecays->At(1))->SetMass(m2);
553 }
554 ///////////////////////////////////////////////////////////////////////////
555 Int_t AliTrack::GetNdecay()
556 {
557 // Provide the number of decay produced tracks
558  return fNdec;
559 }
560 ///////////////////////////////////////////////////////////////////////////
561 AliTrack* AliTrack::GetDecayTrack(Int_t j)
562 {
563 // Provide decay produced track number j
564 // Note : j=1 denotes the first decay track
565  if (!fDecays)
566  {
567   cout << " *AliTrack::GetDecayTrack* No tracks present." << endl;
568   return 0;
569  }
570  else
571  {
572   if ((j >= 1) && (j <= fNdec))
573   {
574    return (AliTrack*)fDecays->At(j-1);
575   }
576   else
577   {
578    cout << " *AliTrack* decay track number : " << j << " out of range."
579         << " Ndec = " << fNdec << endl;
580    return 0;  
581   }
582  }
583 }
584 ///////////////////////////////////////////////////////////////////////////
585 void AliTrack::AddSignal(AliSignal& s)
586 {
587 // Relate an AliSignal object to this track.
588  if (!fSignals) fSignals=new TObjArray();
589  fNsig++;
590  fSignals->Add(&s);
591 }
592 ///////////////////////////////////////////////////////////////////////////
593 void AliTrack::RemoveSignal(AliSignal& s)
594 {
595 // Remove related AliSignal object to this track.
596  if (fSignals)
597  {
598   AliSignal* test=(AliSignal*)fSignals->Remove(&s);
599   if (test)
600   {
601    fNsig--;
602    fSignals->Compress();
603   }
604  }
605 }
606 ///////////////////////////////////////////////////////////////////////////
607 Int_t AliTrack::GetNsignals()
608 {
609 // Provide the number of related AliSignals.
610  return fNsig;
611 }
612 ///////////////////////////////////////////////////////////////////////////
613 AliSignal* AliTrack::GetSignal(Int_t j)
614 {
615 // Provide the related AliSignal number j.
616 // Note : j=1 denotes the first signal.
617  if (!fSignals)
618  {
619   cout << " *AliTrack::GetSignal* No signals present." << endl;
620   return 0;
621  }
622  else
623  {
624   if ((j >= 1) && (j <= fNsig))
625   {
626    return (AliSignal*)fSignals->At(j-1);
627   }
628   else
629   {
630    cout << " *AliTrack* signal number : " << j << " out of range."
631         << " Nsig = " << fNsig << endl;
632    return 0;
633   }
634  }
635 }
636 ///////////////////////////////////////////////////////////////////////////
637 void AliTrack::SetBeginPoint(AliPosition& p)
638 {
639 // Store the position of the track begin-point.
640  if (!fBegin)
641  {
642   fBegin=new AliPositionObj(p);
643  }
644  else
645  {
646   fBegin->Load(p);
647  }
648 }
649 ///////////////////////////////////////////////////////////////////////////
650 AliPosition* AliTrack::GetBeginPoint()
651 {
652 // Provide the position of the track begin-point.
653  return fBegin;
654 }
655 ///////////////////////////////////////////////////////////////////////////
656 void AliTrack::SetEndPoint(AliPosition& p)
657 {
658 // Store the position of the track end-point.
659  if (!fEnd)
660  {
661   fEnd=new AliPositionObj(p);
662  }
663  else
664  {
665   fEnd->Load(p);
666  }
667 }
668 ///////////////////////////////////////////////////////////////////////////
669 AliPosition* AliTrack::GetEndPoint()
670 {
671 // Provide the position of the track end-point.
672  return fEnd;
673 }
674 ///////////////////////////////////////////////////////////////////////////
675 void AliTrack::AddMassHypothesis(Double_t prob,Double_t m,Double_t dm)
676 {
677 // Add a mass hypothesis for this current track.
678 // prob=probalility  m=mass value  dm=error on the mass value.
679 // The default value for the mass error dm is 0.
680  if (!fMasses) fMasses=new TArrayD();
681  if (!fDmasses) fDmasses=new TArrayD();
682  if (!fPmasses) fPmasses=new TArrayD();
683
684  fNmasses++;
685  fMasses->Set(fNmasses);
686  fDmasses->Set(fNmasses);
687  fPmasses->Set(fNmasses);
688
689  fMasses->AddAt(m,fNmasses-1);
690  fDmasses->AddAt(dm,fNmasses-1);
691  fPmasses->AddAt(prob,fNmasses-1);
692 }
693 ///////////////////////////////////////////////////////////////////////////
694 Int_t AliTrack::GetNMassHypotheses()
695 {
696 // Provide the number of mass hypotheses for this track.
697  return fNmasses;
698 }
699 ///////////////////////////////////////////////////////////////////////////
700 Double_t AliTrack::GetMassHypothesis(Int_t j)
701 {
702 // Provide the mass of the jth hypothesis for this track.
703 // Note : the first hypothesis is indicated by j=1.
704 // Default : j=0 ==> Hypothesis with highest probability.
705 // The error on the mass can be obtained by invoking GetResultError()
706 // after invokation of GetMassHypothesis(j).
707
708  Double_t m=0,dm=0,prob=0;
709
710  // Check validity of index j
711  if (j<0 || j>fNmasses)
712  {
713   cout << " *AliTrack::GetMassHypothesis* Invalid index j : " << j
714        << " Number of mass hypotheses : " << fNmasses << endl;
715   fDresult=0;
716   return 0;
717  }
718
719  // Select mass hypothesis with highest probability
720  if (j==0) 
721  {
722   if (fNmasses) 
723   {
724    m=fMasses->At(0);
725    dm=fDmasses->At(0);
726    prob=fPmasses->At(0);
727    for (Int_t i=1; i<fNmasses; i++)
728    {
729     if (fPmasses->At(i)>prob)
730     {
731      m=fMasses->At(i);
732      dm=fDmasses->At(i);
733     }
734    }
735   }
736   fDresult=dm;
737   return m;  
738  }
739
740  // Provide data of requested mass hypothesis
741  m=fMasses->At(j-1);
742  fDresult=fDmasses->At(j-1);
743  return m;
744 }
745 ///////////////////////////////////////////////////////////////////////////
746 Double_t AliTrack::GetMassHypothesisProb(Int_t j)
747 {
748 // Provide the probability of the jth hypothesis for this track.
749 // Note : the first hypothesis is indicated by j=1.
750 // Default : j=0 ==> Hypothesis with highest probability.
751
752  Double_t prob=0;
753
754  // Check validity of index j
755  if (j<0 || j>fNmasses)
756  {
757   cout << " *AliTrack::GetMassHypothesisProb* Invalid index j : " << j
758        << " Number of mass hypotheses : " << fNmasses << endl;
759   return 0;
760  }
761
762  // Select mass hypothesis with highest probability
763  if (j==0) 
764  {
765   if (fNmasses) 
766   {
767    prob=fPmasses->At(0);
768    for (Int_t i=1; i<fNmasses; i++)
769    {
770     if (fPmasses->At(i)>prob) prob=fPmasses->At(i);
771    }
772   }
773   return prob;  
774  }
775
776  // Provide probability of requested mass hypothesis
777  prob=fPmasses->At(j-1);
778  return prob;
779 }
780 ///////////////////////////////////////////////////////////////////////////
781 void AliTrack::SetMass()
782 {
783 // Set the mass and error to the value of the hypothesis with highest prob.
784
785  Double_t m=0,dm=0,prob=0;
786
787  // Select mass hypothesis with highest probability
788  if (fNmasses) 
789  {
790   m=fMasses->At(0);
791   dm=fDmasses->At(0);
792   prob=fPmasses->At(0);
793   for (Int_t i=1; i<fNmasses; i++)
794   {
795    if (fPmasses->At(i)>prob)
796    {
797     m=fMasses->At(i);
798     dm=fDmasses->At(i);
799    }
800   }
801   SetMass(m,dm);
802  }
803  else
804  {
805   cout << " *AliTrack::SetMass()* No hypothesis present => No action." << endl;
806  }
807 }
808 ///////////////////////////////////////////////////////////////////////////
809 void AliTrack::RemoveMassHypothesis(Int_t j)
810 {
811 // Remove the jth mass hypothesis for this track.
812 // Note : the first hypothesis is indicated by j=1.
813
814  if (j<=0 || j>fNmasses) // Check validity of index j
815  {
816   cout << " *AliTrack::RemoveMassHypothesis* Invalid index j : " << j
817        << " Number of mass hypotheses : " << fNmasses << endl;
818  }
819  else
820  {
821   if (j != fNmasses)
822   {
823    fMasses->AddAt(fMasses->At(fNmasses-1),j-1);
824    fDmasses->AddAt(fDmasses->At(fNmasses-1),j-1);
825    fPmasses->AddAt(fPmasses->At(fNmasses-1),j-1);
826   }
827   fMasses->AddAt(0,fNmasses-1);
828   fDmasses->AddAt(0,fNmasses-1);
829   fPmasses->AddAt(0,fNmasses-1);
830   fNmasses--;
831   fMasses->Set(fNmasses);
832   fDmasses->Set(fNmasses);
833   fPmasses->Set(fNmasses);
834  }
835 }
836 ///////////////////////////////////////////////////////////////////////////
837 Double_t AliTrack::GetPt()
838 {
839 // Provide trans. momentum value w.r.t. z-axis.
840 // The error on the value can be obtained by GetResultError()
841 // after invokation of GetPt().
842  Ali3Vector v;
843  v=GetVecTrans();
844  Double_t norm=v.GetNorm();
845  fDresult=v.GetResultError();
846
847  return norm;
848 }
849 ///////////////////////////////////////////////////////////////////////////
850 Double_t AliTrack::GetPl()
851 {
852 // Provide long. momentum value w.r.t. z-axis.
853 // Note : the returned value can also be negative.
854 // The error on the value can be obtained by GetResultError()
855 // after invokation of GetPl().
856  Ali3Vector v;
857  v=GetVecLong();
858
859  Double_t pl=v.GetNorm();
860  fDresult=v.GetResultError();
861
862  Double_t a[3];
863  v.GetVector(a,"sph");
864  if (cos(a[1])<0) pl=-pl;
865
866  return pl;
867 }
868 ///////////////////////////////////////////////////////////////////////////
869 Double_t AliTrack::GetEt()
870 {
871 // Provide trans. energy value w.r.t. z-axis.
872 // The error on the value can be obtained by GetResultError()
873 // after invokation of GetEt().
874  Double_t et=GetScaTrans();
875
876  return et;
877 }
878 ///////////////////////////////////////////////////////////////////////////
879 Double_t AliTrack::GetEl()
880 {
881 // Provide long. energy value w.r.t. z-axis.
882 // Note : the returned value can also be negative.
883 // The error on the value can be obtained by GetResultError()
884 // after invokation of GetEl().
885  Double_t el=GetScaLong();
886
887  return el;
888 }
889 ///////////////////////////////////////////////////////////////////////////
890 Double_t AliTrack::GetMt()
891 {
892 // Provide transverse mass value w.r.t. z-axis.
893 // The error on the value can be obtained by GetResultError()
894 // after invokation of GetMt().
895  Double_t pt=GetPt();
896  Double_t dpt=GetResultError();
897  Double_t m=GetMass();
898  Double_t dm=GetResultError();
899
900  Double_t mt=sqrt(pt*pt+m*m);
901  Double_t dmt2=0;
902  if (mt) dmt2=(pow((pt*dpt),2)+pow((m*dm),2))/(mt*mt);
903
904  fDresult=sqrt(dmt2);
905  return mt;
906 }
907 ///////////////////////////////////////////////////////////////////////////
908 Double_t AliTrack::GetMt(Int_t j)
909 {
910 // Provide transverse mass value w.r.t. z-axis and jth mass hypothesis.
911 // Note : the first hypothesis is indicated by j=1.
912 //        j=0 ==> Hypothesis with highest probability.
913 // The error on the value can be obtained by GetResultError()
914 // after invokation of GetMt(j).
915  Double_t pt=GetPt();
916  Double_t dpt=GetResultError();
917  Double_t m=GetMassHypothesis(j);
918  Double_t dm=GetResultError();
919
920  Double_t mt=sqrt(pt*pt+m*m);
921  Double_t dmt2=0;
922  if (mt) dmt2=(pow((pt*dpt),2)+pow((m*dm),2))/(mt*mt);
923
924  fDresult=sqrt(dmt2);
925  return mt;
926 }
927 ///////////////////////////////////////////////////////////////////////////
928 Double_t AliTrack::GetRapidity()
929 {
930 // Provide rapidity value w.r.t. z-axis.
931 // The error on the value can be obtained by GetResultError()
932 // after invokation of GetRapidity().
933 // Note : Also GetPseudoRapidity() is available since this class is
934 //        derived from Ali4Vector.
935  Double_t e=GetEnergy();
936  Double_t de=GetResultError();
937  Double_t pl=GetPl();
938  Double_t dpl=GetResultError();
939  Double_t sum=e+pl;
940  Double_t dif=e-pl;
941
942  Double_t y=9999,dy2=0;
943  if (sum && dif) y=0.5*log(sum/dif);
944
945  if (sum*dif) dy2=(1./(sum*dif))*(pow((pl*de),2)+pow((e*dpl),2));
946
947  fDresult=sqrt(dy2);
948  return y;
949 }
950 ///////////////////////////////////////////////////////////////////////////
951 void AliTrack::SetImpactPoint(AliPosition& p,TString q)
952 {
953 // Store the position of the impact-point in the plane "q=0".
954 // Here q denotes one of the axes X, Y or Z.
955 // Note : The character to denote the axis may be entered in lower or
956 //        in uppercase.
957  Int_t axis=0;
958  if (q=="x" || q=="X") axis=1;
959  if (q=="y" || q=="Y") axis=2;
960  if (q=="z" || q=="Z") axis=3;
961
962  switch (axis)
963  {
964   case 1: // Impact-point in the plane X=0
965    if (!fImpactYZ)
966    {
967     fImpactYZ=new AliPositionObj(p);
968    }
969    else
970    {
971     fImpactYZ->Load(p);
972    }
973    break;
974
975   case 2: // Impact-point in the plane Y=0
976    if (!fImpactXZ)
977    {
978     fImpactXZ=new AliPositionObj(p);
979    }
980    else
981    {
982     fImpactXZ->Load(p);
983    }
984    break;
985
986   case 3: // Impact-point in the plane Z=0
987    if (!fImpactXY)
988    {
989     fImpactXY=new AliPositionObj(p);
990    }
991    else
992    {
993     fImpactXY->Load(p);
994    }
995    break;
996
997   default: // Unsupported axis
998    cout << "*AliTrack::SetImpactPoint* Unsupported axis : " << q << endl
999         << " Possible axes are 'X', 'Y' and 'Z'." << endl; 
1000    break;
1001  }
1002 }
1003 ///////////////////////////////////////////////////////////////////////////
1004 AliPosition* AliTrack::GetImpactPoint(TString q)
1005 {
1006 // Provide the position of the impact-point in the plane "q=0".
1007 // Here q denotes one of the axes X, Y or Z.
1008 // Note : The character to denote the axis may be entered in lower or
1009 //        in uppercase.
1010  Int_t axis=0;
1011  if (q=="x" || q=="X") axis=1;
1012  if (q=="y" || q=="Y") axis=2;
1013  if (q=="z" || q=="Z") axis=3;
1014
1015  switch (axis)
1016  {
1017   case 1: // Impact-point in the plane X=0
1018    return fImpactYZ;
1019
1020   case 2: // Impact-point in the plane Y=0
1021    return fImpactXZ;
1022
1023   case 3: // Impact-point in the plane Z=0
1024    return fImpactXY;
1025
1026   default: // Unsupported axis
1027    cout << "*AliTrack::GetImpactPoint* Unsupported axis : " << q << endl
1028         << " Possible axes are 'X', 'Y' and 'Z'." << endl; 
1029    return 0;
1030  }
1031 }
1032 ///////////////////////////////////////////////////////////////////////////
1033 void AliTrack::SetId(Int_t id)
1034 {
1035 // Set a user defined unique identifier for this track.
1036  fUserId=id;
1037 }
1038 ///////////////////////////////////////////////////////////////////////////
1039 Int_t AliTrack::GetId()
1040 {
1041 // Provide the user defined unique identifier of this track.
1042  return fUserId;
1043 }
1044 ///////////////////////////////////////////////////////////////////////////
1045 void AliTrack::SetClosestPoint(AliPosition& p)
1046 {
1047 // Set position p as the point of closest approach w.r.t. some reference
1048  if (!fClosest)
1049  {
1050   fClosest=new AliPositionObj(p);
1051  }
1052  else
1053  {
1054   fClosest->Load(p);
1055  }
1056 }
1057 ///////////////////////////////////////////////////////////////////////////
1058 AliPosition* AliTrack::GetClosestPoint()
1059 {
1060 // Provide the point of closest approach w.r.t. some reference
1061  return fClosest;
1062 }
1063 ///////////////////////////////////////////////////////////////////////////
1064 void AliTrack::SetChi2(Float_t chi2)
1065 {
1066 // Set the chi-squared value of the track fit.
1067  if (chi2<0)
1068  {
1069   cout << " *AliTrack::SetChi2* Invalid chi2 value : " << chi2 << endl;
1070  }
1071  else
1072  {
1073   fChi2=chi2;
1074  }
1075 }
1076 ///////////////////////////////////////////////////////////////////////////
1077 void AliTrack::SetNdf(Int_t ndf)
1078 {
1079 // Set the number of degrees of freedom for the track fit.
1080  if (ndf<0)
1081  {
1082   cout << " *AliTrack::SetNdf* Invalid ndf value : " << ndf << endl;
1083  }
1084  else
1085  {
1086   fNdf=ndf;
1087  }
1088 }
1089 ///////////////////////////////////////////////////////////////////////////
1090 Float_t AliTrack::GetChi2()
1091 {
1092 // Provide the chi-squared value of the track fit.
1093  return fChi2;
1094 }
1095 ///////////////////////////////////////////////////////////////////////////
1096 Int_t AliTrack::GetNdf()
1097 {
1098 // Provide the number of degrees of freedom for the track fit.
1099  return fNdf;
1100 }
1101 ///////////////////////////////////////////////////////////////////////////
1102 void AliTrack::SetParticleCode(Int_t code)
1103 {
1104 // Set the user defined particle id code (e.g. the PDF convention).
1105  fCode=code;
1106 }
1107 ///////////////////////////////////////////////////////////////////////////
1108 Int_t AliTrack::GetParticleCode()
1109 {
1110 // Provide the user defined particle id code.
1111  return fCode;
1112 }
1113 ///////////////////////////////////////////////////////////////////////////
1114 void AliTrack::SetParentTrack(AliTrack* t)
1115 {
1116 // Set pointer to the parent track.
1117  fParent=t;
1118 }
1119 ///////////////////////////////////////////////////////////////////////////
1120 AliTrack* AliTrack::GetParentTrack()
1121 {
1122 // Provide pointer to the parent track.
1123  return fParent;
1124 }
1125 ///////////////////////////////////////////////////////////////////////////