]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - EVE/alice-macros/esd_tracks.C
From Pawel Debski.
[u/mrichter/AliRoot.git] / EVE / alice-macros / esd_tracks.C
index 1c6598b1e3849490dc05b77c376aa7da3c0355de..c7aa76cc1d06f6bf12c74a8d9ddfa9d4b04b5db0 100644 (file)
@@ -7,9 +7,65 @@
  * full copyright notice.                                                 *
  **************************************************************************/
 
+#if !defined(__CINT__) || defined(__MAKECINT__)
+
+#include <TString.h>
+#include <TSystem.h>
+#include <TROOT.h>
+#include <TMath.h>
+#include <TGListTree.h>
+
+#include <TEveVSDStructs.h>
+#include <TEveManager.h>
+#include <TEveTrackPropagator.h>
+
+#include <EveBase/AliEveTrack.h>
+#include <EveBase/AliEveTrackCounter.h>
+#include <EveBase/AliEveMagField.h>
+#include <EveBase/AliEveEventManager.h>
+
+#include <AliESDEvent.h>
+#include <AliESDtrackCuts.h>
+#include <AliESDtrack.h>
+#include <AliExternalTrackParam.h>
+
+#include <AliPWG0Helper.h>
+
+#endif
+
+
 // Use inner-tpc track params when its refit failed.
 Bool_t g_esd_tracks_use_ip_on_failed_its_refit = kFALSE;
 
+// Use magnetic-field as retrieved from GRP.
+Bool_t g_esd_tracks_true_field = kTRUE;
+
+// Use Runge-Kutta track stepper.
+Bool_t g_esd_tracks_rk_stepper = kFALSE;
+
+
+//==============================================================================
+
+void esd_track_propagator_setup(TEveTrackPropagator* trkProp,
+                               Float_t magF, Float_t maxR)
+{
+  if (g_esd_tracks_true_field)
+  {
+    trkProp->SetMagFieldObj(new AliEveMagField);
+  }
+  else
+  {
+    trkProp->SetMagField(magF);
+  }
+  if (g_esd_tracks_rk_stepper)
+  {
+    trkProp->SetStepper(TEveTrackPropagator::kRungeKutta);
+  }
+  trkProp->SetMaxR(maxR);
+}
+
+//==============================================================================
+
 TString esd_track_title(AliESDtrack* t)
 {
   TString s;
@@ -21,14 +77,23 @@ TString esd_track_title(AliESDtrack* t)
   Double_t p[3], v[3];
   t->GetXYZ(v);
   t->GetPxPyPz(p);
+  Double_t pt    = t->Pt();
+  Double_t ptsig = TMath::Sqrt(t->GetSigma1Pt2());
+  Double_t ptsq  = pt*pt;
+  Double_t ptm   = pt / (1.0 + pt*ptsig);
+  Double_t ptM   = pt / (1.0 - pt*ptsig);
 
   s = Form("Index=%s, Label=%s\nChg=%d, Pdg=%d\n"
-          "pT=%.3f, pZ=%.3f\nV=(%.3f, %.3f, %.3f)\n",
+          "pT = %.3f + %.3f - %.3f [%.3f]\n"
+           "P  = (%.3f, %.3f, %.3f)\n"
+           "V  = (%.3f, %.3f, %.3f)\n",
           idx.Data(), lbl.Data(), t->Charge(), 0,
-          TMath::Sqrt(p[0]*p[0] + p[1]*p[1]), p[2], v[0], v[1], v[2]);
+          pt, ptM - pt, pt - ptm, ptsig*ptsq,
+           p[0], p[1], p[2],
+           v[0], v[1], v[2]);
 
   Int_t   o;
-  s += "Det(in,out,refit,pid):\n";
+  s += "Det (in,out,refit,pid):\n";
   o  = AliESDtrack::kITSin;
   s += Form("ITS (%d,%d,%d,%d)  ",  t->IsOn(o), t->IsOn(o<<1), t->IsOn(o<<2), t->IsOn(o<<3));
   o  = AliESDtrack::kTPCin;
@@ -51,7 +116,9 @@ TString esd_track_title(AliESDtrack* t)
   return s;
 }
 
-void esd_track_add_param(AliEveTrack* track, AliExternalTrackParam* tp)
+//==============================================================================
+
+void esd_track_add_param(AliEveTrack* track, const AliExternalTrackParam* tp)
 {
   // Add additional track parameters as a path-mark to track.
 
@@ -68,14 +135,17 @@ void esd_track_add_param(AliEveTrack* track, AliExternalTrackParam* tp)
   track->AddPathMark(pm);
 }
 
-TEveTrack* esd_make_track(AliESDtrack *at, TEveTrackList* cont)
+//==============================================================================
+
+AliEveTrack* esd_make_track(AliESDtrack *at, TEveTrackList* cont)
 {
   // Make a standard track representation and put it into given container.
 
   // Choose which parameters to use a track's starting point.
   // If gkFixFailedITSExtr is TRUE (FALSE by default) and
   // if ITS refit failed, take track parameters at inner TPC radius.
-  AliExternalTrackParam* tp = at;
+
+  const AliExternalTrackParam* tp = at;
 
   Bool_t innerTaken = kFALSE;
   if ( ! at->IsOn(AliESDtrack::kITSrefit) && g_esd_tracks_use_ip_on_failed_its_refit)
@@ -84,21 +154,7 @@ TEveTrack* esd_make_track(AliESDtrack *at, TEveTrackList* cont)
     innerTaken = kTRUE;
   }
 
-  Double_t pbuf[3], vbuf[3];
-  tp->GetXYZ(vbuf);
-  tp->GetPxPyPz(pbuf);
-
-  TEveRecTrack rt;
-  rt.fLabel  = at->GetLabel();
-  rt.fIndex  = at->GetID();
-  rt.fStatus = (Int_t) at->GetStatus();
-  rt.fSign   = at->GetSign();
-  rt.fV.Set(vbuf);
-  rt.fP.Set(pbuf);
-  Double_t ep = at->GetP(), mc = at->GetMass();
-  rt.fBeta = ep/TMath::Sqrt(ep*ep + mc*mc);
-
-  AliEveTrack* track = new AliEveTrack(&rt, cont->GetPropagator());
+  AliEveTrack* track = new AliEveTrack(at, cont->GetPropagator());
   track->SetAttLineAttMarker(cont);
   track->SetName(Form("AliEveTrack %d", at->GetID()));
   track->SetElementTitle(esd_track_title(at));
@@ -117,15 +173,20 @@ TEveTrack* esd_make_track(AliESDtrack *at, TEveTrackList* cont)
   return track;
 }
 
+
+//==============================================================================
+// esd_tracks()
+//==============================================================================
+
 TEveTrackList* esd_tracks()
 {
   AliESDEvent* esd = AliEveEventManager::AssertESD();
 
   TEveTrackList* cont = new TEveTrackList("ESD Tracks");
   cont->SetMainColor(6);
-  TEveTrackPropagator* trkProp = cont->GetPropagator();
-  trkProp->SetMagField(0.1*esd->GetMagneticField());
-  trkProp->SetMaxR    (520);
+
+  esd_track_propagator_setup(cont->GetPropagator(),
+                            0.1*esd->GetMagneticField(), 520);
 
   gEve->AddElement(cont);
 
@@ -145,11 +206,12 @@ TEveTrackList* esd_tracks()
   return cont;
 }
 
-TEveElementList* esd_tracks_MI()
+TEveTrackList* esd_tracks_MI()
 {
   AliESDEvent* esd = AliEveEventManager::AssertESD();
 
-  TEveElementList* cont = new TEveElementList("ESD Tracks MI");
+  TEveTrackList* cont = new TEveTrackList("ESD Tracks MI");
+  cont->SetLineColor(5);
   gEve->AddElement(cont);
 
   Int_t count = 0;
@@ -157,13 +219,15 @@ TEveElementList* esd_tracks_MI()
   {
     ++count;
     AliESDtrack* at = esd->GetTrack(n);
-    TEveLine* l = new TEveLine; 
-    l->SetLineColor(5);
-    at->FillPolymarker(l, esd->GetMagneticField(), 0, 250, 5);
+    AliEveTrack* l = new AliEveTrack(at, cont->GetPropagator());
     l->SetName(Form("ESDTrackMI %d", at->GetID()));
     l->SetElementTitle(esd_track_title(at));
+    l->SetAttLineAttMarker(cont);
     l->SetSourceObject(at);
 
+    at->FillPolymarker(l, esd->GetMagneticField(), 0, 250, 5);
+
+    l->SetLockPoints(kTRUE);
     cont->AddElement(l);
   }
   cont->SetTitle(Form("N=%d", count));
@@ -173,9 +237,10 @@ TEveElementList* esd_tracks_MI()
   return cont;
 }
 
-/******************************************************************************/
+
+//==============================================================================
 // esd_tracks_from_array()
-/******************************************************************************/
+//==============================================================================
 
 TEveTrackList* esd_tracks_from_array(TCollection* col, AliESDEvent* esd=0)
 {
@@ -186,9 +251,9 @@ TEveTrackList* esd_tracks_from_array(TCollection* col, AliESDEvent* esd=0)
 
   TEveTrackList* cont = new TEveTrackList("ESD Tracks");
   cont->SetMainColor(6);
-  TEveTrackPropagator* trkProp = cont->GetPropagator();
-  trkProp->SetMagField(0.1*esd->GetMagneticField());
-  trkProp->SetMaxR    (520);
+
+  esd_track_propagator_setup(cont->GetPropagator(),
+                            0.1*esd->GetMagneticField(), 520);
 
   gEve->AddElement(cont);
 
@@ -205,7 +270,7 @@ TEveTrackList* esd_tracks_from_array(TCollection* col, AliESDEvent* esd=0)
     }
 
     ++count;
-    AliESDtrack* at = (AliESDtrack*) obj;
+    AliESDtrack* at = reinterpret_cast<AliESDtrack*>(obj);
 
     AliEveTrack* track = esd_make_track(at, cont);
     cont->AddElement(track);
@@ -221,18 +286,18 @@ TEveTrackList* esd_tracks_from_array(TCollection* col, AliESDEvent* esd=0)
 void esd_tracks_alianalcuts_demo()
 {
   AliESDEvent* esd = AliEveEventManager::AssertESD();
-  gSystem->Load("libANALYSIS");
 
-  AliAnalysisTrackCuts atc;
+  AliESDtrackCuts atc;
   atc.SetPtRange(0.1, 5);
   atc.SetRapRange(-1, 1);
 
-  esd_tracks_from_array(atc.GetAcceptedParticles(esd), esd);
+  esd_tracks_from_array(atc.GetAcceptedTracks(esd), esd);
 }
 
-/******************************************************************************/
+
+//==============================================================================
 // esd_tracks_by_category
-/******************************************************************************/
+//==============================================================================
 
 Float_t get_sigma_to_vertex(AliESDtrack* esdTrack)
 {
@@ -276,21 +341,186 @@ Float_t get_sigma_to_vertex(AliESDtrack* esdTrack)
   return d;
 }
 
-TEveElementList* g_esd_tracks_by_category_container = 0;
-
 TEveElementList* esd_tracks_by_category()
 {
-  // Import ESD tracks, separate them into several containers
-  // according to primary-vertex cut and ITS&TPC refit status.
+   // Import ESD tracks, separate them into several containers
+   // according to primary-vertex cut and ITS&TPC refit status.
+
+   AliESDEvent* esd = AliEveEventManager::AssertESD();
+
+   TEveElementList* cont = new TEveElementList("ESD Tracks by category");
+   gEve->AddElement(cont);
+
+   const Int_t   nCont = 9;
+   const Float_t maxR  = 520;
+   const Float_t magF  = 0.1*esd->GetMagneticField();
+
+   TEveTrackList *tl[nCont];
+   Int_t          tc[nCont];
+   Int_t          count = 0;
+
+   tl[0] = new TEveTrackList("Sigma < 3");
+   tc[0] = 0;
+   esd_track_propagator_setup(tl[0]->GetPropagator(), magF, maxR);
+   tl[0]->SetMainColor(3);
+   cont->AddElement(tl[0]);
+
+   tl[1] = new TEveTrackList("3 < Sigma < 5");
+   tc[1] = 0;
+   esd_track_propagator_setup(tl[1]->GetPropagator(), magF, maxR);
+   tl[1]->SetMainColor(7);
+   cont->AddElement(tl[1]);
+
+   tl[2] = new TEveTrackList("5 < Sigma");
+   tc[2] = 0;
+   esd_track_propagator_setup(tl[2]->GetPropagator(), magF, maxR);
+   tl[2]->SetMainColor(46);
+   cont->AddElement(tl[2]);
+
+   tl[3] = new TEveTrackList("no ITS refit; Sigma < 5");
+   tc[3] = 0;
+   esd_track_propagator_setup(tl[3]->GetPropagator(), magF, maxR);
+   tl[3]->SetMainColor(41);
+   cont->AddElement(tl[3]);
+
+   tl[4] = new TEveTrackList("no ITS refit; Sigma > 5");
+   tc[4] = 0;
+   esd_track_propagator_setup(tl[4]->GetPropagator(), magF, maxR);
+   tl[4]->SetMainColor(48);
+   cont->AddElement(tl[4]);
+
+   tl[5] = new TEveTrackList("no TPC refit");
+   tc[5] = 0;
+   esd_track_propagator_setup(tl[5]->GetPropagator(), magF, maxR);
+   tl[5]->SetMainColor(kRed);
+   cont->AddElement(tl[5]);
+
+   tl[6] = new TEveTrackList("ITS ncl>=3 & SPD Inner");
+   tc[6] = 0;
+   esd_track_propagator_setup(tl[6]->GetPropagator(), magF, maxR);
+   tl[6]->SetMainColor(kGreen);
+   cont->AddElement(tl[6]);
+
+   tl[7] = new TEveTrackList("ITS ncl>=3 & b<3 cm");
+   tc[7] = 0;
+   esd_track_propagator_setup(tl[7]->GetPropagator(), magF, maxR);
+   tl[7]->SetMainColor(kMagenta - 9);
+   cont->AddElement(tl[7]);
+
+   tl[8] = new TEveTrackList("ITS others");
+   tc[8] = 0;
+   esd_track_propagator_setup(tl[8]->GetPropagator(), magF, maxR);
+   tl[8]->SetMainColor(kRed + 2);
+   cont->AddElement(tl[8]);
+
+   for (Int_t n = 0; n < esd->GetNumberOfTracks(); ++n)
+   {
+      AliESDtrack* at = esd->GetTrack(n);
+
+      Float_t s  = get_sigma_to_vertex(at);
+      Int_t   ti;
+      if      (s <  3) ti = 0;
+      else if (s <= 5) ti = 1;
+      else             ti = 2;
+
+      Int_t    nclits;
+      Double_t dtobeam;
+
+      if (at->IsOn(AliESDtrack::kITSin) && ! at->IsOn(AliESDtrack::kTPCin))
+      {
+         UChar_t itsclmap = at->GetITSClusterMap();
+         Bool_t  spdinner = (itsclmap & 3) != 0;
+
+         nclits = 0;
+         for (Int_t iter = 0; iter < 6; ++iter)
+            if (itsclmap & (1 << iter)) nclits ++;
+
+         Double_t xyz[3];
+         at->GetXYZ(xyz);
+         dtobeam = TMath::Hypot(xyz[0], xyz[1]);
+
+         if ((nclits >= 3) && (spdinner))
+            ti = 6;
+         else if ((nclits >= 3) && (dtobeam < 3.0))
+            ti = 7;
+         else
+            ti = 8;
+      }
+      else if (at->IsOn(AliESDtrack::kTPCin) && ! at->IsOn(AliESDtrack::kTPCrefit))
+      {
+         ti = 5;
+      }
+      else if ( ! at->IsOn(AliESDtrack::kITSrefit))
+      {
+         ti = (ti == 2) ? 4 : 3;
+      }
+
+      TEveTrackList* tlist = tl[ti];
+      ++tc[ti];
+      ++count;
+
+      AliEveTrack* track = esd_make_track(at, tlist);
+      if (ti<6)
+         track->SetName(Form("ESD Track idx=%d, sigma=%5.3f", at->GetID(), s));
+      else
+         track->SetName(Form("ESD Track idx=%d, dxy=%5.3f cl=%i", at->GetID(), dtobeam, nclits));
+      tlist->AddElement(track);
+   }
+
+   for (Int_t ti = 0; ti < nCont; ++ti)
+   {
+      TEveTrackList* tlist = tl[ti];
+      tlist->SetName(Form("%s [%d]", tlist->GetName(), tlist->NumChildren()));
+      tlist->SetTitle(Form("N tracks=%d", tc[ti]));
+
+      tlist->MakeTracks();
+
+      //    Bool_t good_cont = ti <= 1;
+      Bool_t good_cont = ((ti == 6) || (ti == 7));
+      if (AliEveTrackCounter::IsActive())
+      {
+         AliEveTrackCounter::fgInstance->RegisterTracks(tlist, good_cont);
+      }
+      else
+      {
+         if ( ! good_cont)
+            tlist->SetLineStyle(6);
+      }
+   }
+   cont->SetTitle(Form("N all tracks = %d", count));
+   // ??? The following does not always work:
+   cont->FindListTreeItem(gEve->GetListTree())->SetOpen(kTRUE);
+
+   gEve->Redraw3D();
+
+   return cont;
+}
+
+
+//==============================================================================
+// esd_tracks_by_anal_cuts
+//==============================================================================
+
+AliESDtrackCuts* g_esd_tracks_anal_cuts = 0;
 
+TEveElementList* esd_tracks_by_anal_cuts()
+{
   AliESDEvent* esd = AliEveEventManager::AssertESD();
 
-  TEveElementList* cont = new TEveElementList("ESD Tracks by category");
-  gEve->AddElement(cont);
+  if (g_esd_tracks_anal_cuts == 0)
+  {
+    gSystem->Load("libPWG0base");
+    gROOT->ProcessLine(".L $ALICE_ROOT/PWG0/CreateStandardCuts.C");
+    Int_t mode = AliPWG0Helper::kTPC;
+    if (TMath::Abs(esd->GetMagneticField()) > 0.01)
+      mode |= AliPWG0Helper::kFieldOn;
+    gROOT->ProcessLine(Form("g_esd_tracks_anal_cuts = CreateTrackCuts(%d, kFALSE)", mode));
+  }
 
-  g_esd_tracks_by_category_container = cont;
+  TEveElementList* cont = new TEveElementList("ESD Tracks by Analysis Cuts");
+  gEve->AddElement(cont);
 
-  const Int_t   nCont = 7;
+  const Int_t   nCont = 2;
   const Float_t maxR  = 520;
   const Float_t magF  = 0.1*esd->GetMagneticField();
 
@@ -298,84 +528,31 @@ TEveElementList* esd_tracks_by_category()
   Int_t          tc[nCont];
   Int_t          count = 0;
 
-  tl[0] = new TEveTrackList("Sigma < 3");
+  tl[0] = new TEveTrackList("Passed");
   tc[0] = 0;
-  tl[0]->GetPropagator()->SetMagField(magF);
-  tl[0]->GetPropagator()->SetMaxR    (maxR);
+  esd_track_propagator_setup(tl[0]->GetPropagator(), magF, maxR);
   tl[0]->SetMainColor(3);
   cont->AddElement(tl[0]);
 
-  tl[1] = new TEveTrackList("3 < Sigma < 5");
+  tl[1] = new TEveTrackList("Rejected");
   tc[1] = 0;
-  tl[1]->GetPropagator()->SetMagField(magF);
-  tl[1]->GetPropagator()->SetMaxR    (maxR);
-  tl[1]->SetMainColor(7);
+  esd_track_propagator_setup(tl[1]->GetPropagator(), magF, maxR);
+  tl[1]->SetMainColor(kRed);
   cont->AddElement(tl[1]);
 
-  tl[2] = new TEveTrackList("5 < Sigma");
-  tc[2] = 0;
-  tl[2]->GetPropagator()->SetMagField(magF);
-  tl[2]->GetPropagator()->SetMaxR    (maxR);
-  tl[2]->SetMainColor(46);
-  cont->AddElement(tl[2]);
-
-  tl[3] = new TEveTrackList("no ITS refit; Sigma < 5");
-  tc[3] = 0;
-  tl[3]->GetPropagator()->SetMagField(magF);
-  tl[3]->GetPropagator()->SetMaxR    (maxR);
-  tl[3]->SetMainColor(41);
-  cont->AddElement(tl[3]);
-
-  tl[4] = new TEveTrackList("no ITS refit; Sigma > 5");
-  tc[4] = 0;
-  tl[4]->GetPropagator()->SetMagField(magF);
-  tl[4]->GetPropagator()->SetMaxR    (maxR);
-  tl[4]->SetMainColor(48);
-  cont->AddElement(tl[4]);
-
-  tl[5] = new TEveTrackList("no TPC refit");
-  tc[5] = 0;
-  tl[5]->GetPropagator()->SetMagField(magF);
-  tl[5]->GetPropagator()->SetMaxR    (maxR);
-  tl[5]->SetMainColor(kRed);
-  cont->AddElement(tl[5]);
-
-  tl[6] = new TEveTrackList("ITS stand-alone");
-  tc[6] = 0;
-  tl[6]->GetPropagator()->SetMagField(magF);
-  tl[6]->GetPropagator()->SetMaxR    (maxR);
-  tl[6]->SetMainColor(kMagenta - 9);
-  cont->AddElement(tl[6]);
-
   for (Int_t n = 0; n < esd->GetNumberOfTracks(); ++n)
   {
     AliESDtrack* at = esd->GetTrack(n);
 
     Float_t s  = get_sigma_to_vertex(at);
-    Int_t   ti;
-    if      (s <  3) ti = 0;
-    else if (s <= 5) ti = 1;
-    else             ti = 2;
-
-    if (at->IsOn(AliESDtrack::kITSin) && ! at->IsOn(AliESDtrack::kTPCin))
-    {
-      ti = 6;
-    }
-    else if (at->IsOn(AliESDtrack::kTPCin) && ! at->IsOn(AliESDtrack::kTPCrefit))
-    {
-      ti = 5;
-    }
-    else if ( ! at->IsOn(AliESDtrack::kITSrefit))
-    {
-      ti = (ti == 2) ? 4 : 3;
-    }
+    Int_t   ti = (g_esd_tracks_anal_cuts->AcceptTrack(at)) ? 0 : 1;
 
     TEveTrackList* tlist = tl[ti];
     ++tc[ti];
     ++count;
 
     AliEveTrack* track = esd_make_track(at, tlist);
-    track->SetName(Form("AliEveTrack idx=%d, sigma=%5.3f", at->GetID(), s));
+    track->SetName(Form("ESD Track idx=%d, sigma=%5.3f", at->GetID(), s));
     tlist->AddElement(track);
   }
 
@@ -386,6 +563,17 @@ TEveElementList* esd_tracks_by_category()
     tlist->SetTitle(Form("N tracks=%d", tc[ti]));
 
     tlist->MakeTracks();
+
+    Bool_t good_cont = ti < 1;
+    if (AliEveTrackCounter::IsActive())
+    {
+      AliEveTrackCounter::fgInstance->RegisterTracks(tlist, good_cont);
+    }
+    else
+    {
+      if ( ! good_cont)
+        tlist->SetLineStyle(6);
+    }
   }
   cont->SetTitle(Form("N all tracks = %d", count));
   // ??? The following does not always work: