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 ///////////////////////////////////////////////////////////////////////////