fe4f949582d4896074735900c721839c18d93d1a
[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 /*
17 $Log$
18 Revision 1.2  1999/09/29 09:24:28  fca
19 Introduction of the Copyright and cvs Log
20
21 */
22
23 ///////////////////////////////////////////////////////////////////////////
24 // Class AliTrack
25 // Handling of the attributes of a reconstructed particle track.
26 //
27 // Coding example :
28 // ----------------
29 //
30 // Float_t a[4]={195.,1.2,-0.04,8.5};
31 // Ali4Vector pmu;
32 // pmu.SetVector(a,"car");
33 // AliTrack t1;
34 // t1.Set4Momentum(pmu);
35 //
36 // Float_t b[3]={1.2,-0.04,8.5};
37 // Ali3Vector p;
38 // p.SetVector(b,"car");
39 // AliTrack t2;
40 // t2.Set3Momentum(p);
41 // t2.SetCharge(0);
42 // t2.SetMass(1.115);
43 //
44 // t1.Info();
45 // t2.Info();
46 //
47 // Float_t pi=acos(-1.);
48 // Float_t thcms=0.2*pi; // decay theta angle in cms
49 // Float_t phicms=pi/4.; // decay theta angle in cms
50 // Float_t m1=0.938;
51 // Float_t m2=0.140;
52 // t2.Decay(m1,m2,thcms,phicms); // Track t2 decay : Lambda -> proton + pion
53 //
54 // t2.List();
55 //
56 // Int_t ndec=t2.GetNdecay();
57 // AliTrack* d1=t2.GetDecayTrack(1); // Access to decay track number 1
58 // AliTrack* d2=t2.GetDecayTrack(2); // Access to decay track number 2
59 //
60 // AliSignal s1,s2,s3,s4;
61 //
62 // .... // Code (e.g. detector readout) to fill AliSignal data
63 //
64 // AliTrack trec; // Track which will be reconstructed from signals
65 // trec.AddSignal(s1);
66 // trec.AddSignal(s3);
67 // trec.AddSignal(s4);
68 //
69 // Ali3Vector P;
70 // Float_t Q,M;
71 //
72 // ... // Code which accesses signals from trec and reconstructs
73 //        3-momentum P, charge Q, mass M etc...
74 //
75 // trec.Set3Momentum(P);
76 // trec.SetCharge(Q);
77 // trec.SetMass(M);
78 //
79 // Float_t r1[3]={1.6,-3.8,25.7};
80 // Float_t er1[3]={0.2,0.5,1.8};
81 // Float_t r2[3]={8.6,23.8,-6.7};
82 // Float_t er2[3]={0.93,1.78,0.8};
83 // AliPosition begin,end;
84 // begin.SetPosition(r1,"car");
85 // begin.SetPositionErrors(er1,"car");
86 // end.SetPosition(r2,"car");
87 // end.SetPositionErrors(er2,"car");
88 // trec.SetBeginPoint(begin);
89 // trec.SetEndPoint(end);
90 // 
91 // Note : All quantities are in GeV, GeV/c or GeV/c**2
92 //
93 //--- Author: Nick van Eijndhoven 10-jul-1997 UU-SAP Utrecht
94 //- Modified: NvE 29-oct-1999 UU-SAP Utrecht
95 ///////////////////////////////////////////////////////////////////////////
96
97 #include "AliTrack.h"
98  
99 ClassImp(AliTrack) // Class implementation to enable ROOT I/O
100  
101 AliTrack::AliTrack()
102 {
103 // Default constructor
104 // All variables initialised to 0
105  fDecays=0;
106  fSignals=0;
107  Reset();
108 }
109 ///////////////////////////////////////////////////////////////////////////
110 AliTrack::~AliTrack()
111 {
112 // Destructor to delete memory allocated for decay tracks array
113  if (fDecays)
114  {
115   fDecays->Delete();
116   delete fDecays;
117   fDecays=0;
118  }
119  if (fSignals)
120  {
121   fSignals->Clear();
122   delete fSignals;
123   fSignals=0;
124  }
125 }
126 ///////////////////////////////////////////////////////////////////////////
127 void AliTrack::Reset()
128 {
129 // Reset all variables to 0
130  fQ=0;
131  fNdec=0;
132  fNsig=0;
133  Double_t a[4]={0,0,0,0};
134  SetVector(a,"sph");
135  if (fDecays)
136  {
137   fDecays->Delete();
138   delete fDecays;
139   fDecays=0;
140  }
141  if (fSignals)
142  {
143   fSignals->Clear();
144   delete fSignals;
145   fSignals=0;
146  }
147  Double_t b[3]={0,0,0};
148  fBegin.SetPosition(b,"sph");
149  fEnd.SetPosition(b,"sph");
150 }
151 ///////////////////////////////////////////////////////////////////////////
152 void AliTrack::Set3Momentum(Ali3Vector& p)
153 {
154 // Set the track parameters according to the 3-momentum p
155  Set3Vector(p);
156 }
157 ///////////////////////////////////////////////////////////////////////////
158 void AliTrack::Set4Momentum(Ali4Vector& p)
159 {
160 // Set the track parameters according to the 4-momentum p
161  Double_t E=p.GetScalar();
162  Double_t dE=p.GetResultError();
163  Ali3Vector pv=p.Get3Vector();
164  SetVector(E,pv);
165  SetScalarError(dE);
166 }
167 ///////////////////////////////////////////////////////////////////////////
168 void AliTrack::SetMass(Double_t m,Double_t dm)
169 {
170 // Set the particle mass
171 // The default value for the error dm is 0.
172  Double_t inv=pow(m,2);
173  Double_t dinv=fabs(2.*m*dm);
174  SetInvariant(inv,dinv);
175 }
176 ///////////////////////////////////////////////////////////////////////////
177 void AliTrack::SetCharge(Float_t q)
178 {
179 // Set the particle charge
180  fQ=q;
181 }
182 ///////////////////////////////////////////////////////////////////////////
183 void AliTrack::Info(TString f)
184 {
185 // Provide track information within the coordinate frame f
186  Double_t m=GetMass();
187  Double_t dm=GetResultError();
188  cout << " *AliTrack::Info* Mass : " << m
189       << " error : " << dm << " Charge : " << fQ
190       << " Momentum : " << GetMomentum() << " Ntracks : " << fNdec
191       << " Nsignals : " << fNsig << endl;
192  Ali4Vector::Info(f); 
193
194 ///////////////////////////////////////////////////////////////////////////
195 void AliTrack::List(TString f)
196 {
197 // Provide current track and decay level 1 information within coordinate frame f
198
199  Info(f); // Information of the current track
200
201  // Decay products of this track
202  AliTrack* td; 
203  for (Int_t id=1; id<=fNdec; id++)
204  {
205   td=GetDecayTrack(id);
206   if (td)
207   {
208    cout << "  ---Level 1 sec. track no. " << id << endl;
209    td->Info(f); 
210   }
211   else
212   {
213    cout << " *AliTrack::List* Error : No decay track present." << endl; 
214   }
215  }
216
217 ///////////////////////////////////////////////////////////////////////////
218 void AliTrack::ListAll(TString f)
219 {
220 // Provide complete track and decay information within the coordinate frame f
221
222  Info(f); // Information of the current track
223  cout << " Begin-point :"; fBegin.Info(f);
224  cout << " End-point   :"; fEnd.Info(f);
225  for (Int_t is=1; is<=GetNsignals(); is++)
226  {
227   ((AliSignal*)GetSignal(is))->Info(f);
228  }
229
230  AliTrack* t=this;
231  Dump(t,1,f); // Information of all decay products
232 }
233 //////////////////////////////////////////////////////////////////////////
234 void AliTrack::Dump(AliTrack* t,Int_t n,TString f)
235 {
236 // Recursively provide the info of all decay levels of this track
237  AliTrack* td; 
238  for (Int_t id=1; id<=t->GetNdecay(); id++)
239  {
240   td=t->GetDecayTrack(id);
241   if (td)
242   {
243    cout << "  ---Level " << n << " sec. track no. " << id << endl;
244    td->Info(f); 
245    for (Int_t is=1; is<=td->GetNsignals(); is++)
246    {
247     ((AliSignal*)td->GetSignal(is))->Info(f);
248    }
249
250    // Go for next decay level of this decay track recursively
251    Dump(td,n+1,f);
252   }
253   else
254   {
255    cout << " *AliTrack::Dump* Error : No decay track present." << endl; 
256   }
257  }
258
259 //////////////////////////////////////////////////////////////////////////
260 Double_t AliTrack::GetMomentum()
261 {
262 // Provide the value of the track 3-momentum.
263 // The error can be obtained by invoking GetResultError() after
264 // invokation of GetMomentum().
265
266 // Ali3Vector p=Get3Vector();
267 // return sqrt(p.Dot(p));
268  Double_t norm=fV.GetNorm();
269  return norm;
270 }
271 ///////////////////////////////////////////////////////////////////////////
272 Ali3Vector AliTrack::Get3Momentum()
273 {
274 // Provide the track 3-momentum
275  return (Ali3Vector)Get3Vector();
276 }
277 ///////////////////////////////////////////////////////////////////////////
278 Double_t AliTrack::GetMass()
279 {
280 // Provide the particle mass.
281 // The error can be obtained by invoking GetResultError() after
282 // invokation of GetMass().
283  Double_t inv=GetInvariant();
284  Double_t dinv=GetResultError();
285  Double_t dm=0;
286  if (inv >= 0)
287  {
288  Double_t m=sqrt(inv);
289  if (m) dm=dinv/(2.*m);
290  fDresult=dm;
291  return m;
292  }
293  else
294  {
295   cout << "*AliTrack::GetMass* Unphysical situation m**2 = " << inv << endl;
296   cout << " Value 0 will be returned." << endl;
297   fDresult=dm;
298   return 0;
299  }
300 }
301 ///////////////////////////////////////////////////////////////////////////
302 Float_t AliTrack::GetCharge()
303 {
304 // Provide the particle charge
305  return fQ;
306 }
307 ///////////////////////////////////////////////////////////////////////////
308 Double_t AliTrack::GetEnergy()
309 {
310 // Provide the particle's energy.
311 // The error can be obtained by invoking GetResultError() after
312 // invokation of GetEnergy().
313  Double_t E=GetScalar();
314  if (E>0)
315  {
316   return E;
317  }
318  else
319  {
320   cout << "*AliTrack::GetEnergy* Unphysical situation E = " << E << endl;
321   cout << " Value 0 will be returned." << endl;
322   return 0;
323  }
324 }
325 ///////////////////////////////////////////////////////////////////////////
326 void AliTrack::Decay(Double_t m1,Double_t m2,Double_t thcms,Double_t phicms)
327 {
328 // Perform 2-body decay of current track
329 // m1     : mass of decay product 1
330 // m2     : mass of decay product 2
331 // thcms  : cms theta decay angle (in rad.) of m1
332 // phicms : cms phi decay angle (in rad.) of m1
333  
334  fNdec=2; // it's a 2-body decay
335
336  Double_t M=GetMass();
337  
338 // Compute the 4-momenta of the decay products in the cms
339 // Note : p2=p1=pnorm for a 2-body decay
340  Double_t e1=0;
341  if (M) e1=((M*M)+(m1*m1)-(m2*m2))/(2.*M);
342  Double_t e2=0;
343  if (M) e2=((M*M)+(m2*m2)-(m1*m1))/(2.*M);
344  Double_t pnorm=(e1*e1)-(m1*m1);
345  if (pnorm>0.)
346  {
347   pnorm=sqrt(pnorm);
348  }
349  else
350  {
351   pnorm=0;
352  }
353  
354  Double_t a[3];
355  a[0]=pnorm;
356  a[1]=thcms;
357  a[2]=phicms;
358  Ali3Vector p;
359  p.SetVector(a,"sph");
360
361  Ali4Vector pprim1;
362  pprim1.SetVector(e1,p);
363  pprim1.SetInvariant(m1*m1);
364
365  Ali4Vector pprim2;
366  p*=-1;
367  pprim2.SetVector(e2,p);
368  pprim2.SetInvariant(m2*m2);
369
370  // Determine boost parameters from the parent particle
371  Double_t E=GetEnergy();
372  p=Get3Vector();
373  Ali4Vector pmu;
374  pmu.SetVector(E,p);
375
376  AliBoost q;
377  q.Set4Momentum(pmu);
378  
379  Ali4Vector p1=q.Inverse(pprim1); // Boost decay product 1
380  Ali4Vector p2=q.Inverse(pprim2); // Boost decay product 2
381  
382  // Enter the boosted data into the decay tracks array
383  if (fDecays)
384  {
385   fDecays->Delete();
386   delete fDecays;
387  }
388  fDecays=new TObjArray();
389
390  fDecays->Add(new AliTrack);
391  ((AliTrack*)fDecays->At(0))->Set4Momentum(p1);
392  ((AliTrack*)fDecays->At(0))->SetMass(m1);
393  fDecays->Add(new AliTrack);
394  ((AliTrack*)fDecays->At(1))->Set4Momentum(p2);
395  ((AliTrack*)fDecays->At(1))->SetMass(m2);
396 }
397 ///////////////////////////////////////////////////////////////////////////
398 Int_t AliTrack::GetNdecay()
399 {
400 // Provide the number of decay produced tracks
401  return fNdec;
402 }
403 ///////////////////////////////////////////////////////////////////////////
404 AliTrack* AliTrack::GetDecayTrack(Int_t j)
405 {
406 // Provide decay produced track number j
407 // Note : j=1 denotes the first decay track
408  if ((j >= 1) && (j <= fNdec))
409  {
410   return (AliTrack*)fDecays->At(j-1);
411  }
412  else
413  {
414   cout << " *AliTrack* decay track number : " << j << " out of range." << endl;
415   cout << " -- Decay track number 1 (if any) returned." << endl;
416   return (AliTrack*)fDecays->At(0);
417  }
418 }
419 ///////////////////////////////////////////////////////////////////////////
420 void AliTrack::AddSignal(AliSignal& s)
421 {
422 // Relate an AliSignal object to this track.
423  if (!fSignals) fSignals=new TObjArray();
424  fNsig++;
425  fSignals->Add(&s);
426 }
427 ///////////////////////////////////////////////////////////////////////////
428 void AliTrack::RemoveSignal(AliSignal& s)
429 {
430 // Remove related AliSignal object to this track.
431  if (fSignals)
432  {
433   AliSignal* test=(AliSignal*)fSignals->Remove(&s);
434   if (test)
435   {
436    fNsig--;
437    fSignals->Compress();
438   }
439  }
440 }
441 ///////////////////////////////////////////////////////////////////////////
442 Int_t AliTrack::GetNsignals()
443 {
444 // Provide the number of related AliSignals.
445  return fNsig;
446 }
447 ///////////////////////////////////////////////////////////////////////////
448 AliSignal* AliTrack::GetSignal(Int_t j)
449 {
450 // Provide the related AliSignal number j.
451 // Note : j=1 denotes the first signal.
452  if ((j >= 1) && (j <= fNsig))
453  {
454   return (AliSignal*)fSignals->At(j-1);
455  }
456  else
457  {
458   cout << " *AliTrack* signal number : " << j << " out of range." << endl;
459   cout << " -- Signal number 1 (if any) returned." << endl;
460   return (AliSignal*)fDecays->At(0);
461  }
462 }
463 ///////////////////////////////////////////////////////////////////////////
464 void AliTrack::SetBeginPoint(AliPosition p)
465 {
466 // Store the position of the track begin-point.
467  fBegin=p;
468 }
469 ///////////////////////////////////////////////////////////////////////////
470 AliPosition AliTrack::GetBeginPoint()
471 {
472 // Provide the position of the track begin-point.
473  return fBegin;
474 }
475 ///////////////////////////////////////////////////////////////////////////
476 void AliTrack::SetEndPoint(AliPosition p)
477 {
478 // Store the position of the track end-point.
479  fEnd=p;
480 }
481 ///////////////////////////////////////////////////////////////////////////
482 AliPosition AliTrack::GetEndPoint()
483 {
484 // Provide the position of the track end-point.
485  return fEnd;
486 }
487 ///////////////////////////////////////////////////////////////////////////