PWG4/UserTasks/EmcalTasks/AliEmcalEsdTpcTrackTask.cxx

   1 // $Id$
   2 //
   3 // Task to constrain TPC tracks to SPD vertex.
   4 //
   5 //
   6
   7 #include <TClonesArray.h>
   8 #include <TGeoGlobalMagField.h>
   9 #include "AliAnalysisManager.h"
  10 #include "AliESDEvent.h"
  11 #include "AliESDtrackCuts.h"
  12 #include "AliEmcalEsdTpcTrackTask.h"
  13 #include "AliMagF.h"
  14
  15 ClassImp(AliEmcalEsdTpcTrackTask)
  16
  17 //________________________________________________________________________
  18 AliEmcalEsdTpcTrackTask::AliEmcalEsdTpcTrackTask() :
  19   AliAnalysisTaskSE(),
  20   fEsdTrackCuts(0),
  21   fHybridTrackCuts(0),
  22   fTracksName(),
  23   fEsdEv(0),
  24   fTracks(0)
  25 {
  26   // Constructor.
  27 }
  28
  29 //________________________________________________________________________
  30 AliEmcalEsdTpcTrackTask::AliEmcalEsdTpcTrackTask(const char *name) :
  31   AliAnalysisTaskSE(name),
  32   fEsdTrackCuts(0),
  33   fHybridTrackCuts(0),
  34   fTracksName("TpcSpdVertexConstrainedTracks"),
  35   fEsdEv(0),
  36   fTracks(0)
  37 {
  38   // Constructor.
  39   fBranchNames = "ESD:AliESDHeader.,AliESDRun.,SPDVertex.,Tracks";
  40 }
  41
  42 //________________________________________________________________________
  43 AliEmcalEsdTpcTrackTask::~AliEmcalEsdTpcTrackTask()
  44 {
  45   //Destructor
  46
  47   delete fEsdTrackCuts;
  48 }
  49
  50 //________________________________________________________________________
  51 void AliEmcalEsdTpcTrackTask::UserCreateOutputObjects()
  52 {
  53   // Create histograms.
  54
  55   fTracks = new TClonesArray("AliESDtrack");
  56   fTracks->SetName(fTracksName);
  57
  58   if (!fEsdTrackCuts) {
  59     AliInfo("No track cuts given, creating default cuts");
  60     fEsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
  61     fEsdTrackCuts->SetPtRange(0.15,1e3);
  62   }
  63 }
  64
  65 //________________________________________________________________________
  66 void AliEmcalEsdTpcTrackTask::UserExec(Option_t *)
  67 {
  68   // Main loop, called for each event.
  69
  70   fEsdEv = dynamic_cast<AliESDEvent*>(InputEvent());
  71   if (!fEsdEv) {
  72     AliError("Task works only on ESD events, returning");
  73     return;
  74   }
  75
  76   AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
  77   if (!am) {
  78     AliError("Manager zero, returning");
  79     return;
  80   }
  81
  82   // add tracks to event if not yet there
  83   if (!(InputEvent()->FindListObject(fTracksName)))
  84     InputEvent()->AddObject(fTracks);
  85
  86   if (!fHybridTrackCuts) {
  87     am->LoadBranch("AliESDRun.");
  88     am->LoadBranch("AliESDHeader.");
  89     if (!TGeoGlobalMagField::Instance()->GetField()) { // construct field map
  90       fEsdEv->InitMagneticField();
  91     }
  92
  93     am->LoadBranch("SPDVertex.");
  94     const AliESDVertex *vtxSPD = fEsdEv->GetPrimaryVertexSPD();
  95     if (!vtxSPD) {
  96       AliError("No SPD vertex, returning");
  97       return;
  98     }
  99
 100     am->LoadBranch("Tracks");
 101     Int_t ntr = fEsdEv->GetNumberOfTracks();
 102     for (Int_t i=0, ntrnew=0; i<ntr; ++i) {
 103       AliESDtrack *etrack = fEsdEv->GetTrack(i);
 104       if (!etrack)
 105         continue;
 106       if (!fEsdTrackCuts->AcceptTrack(etrack))
 107         continue;
 108
 109       AliESDtrack *ntrack = AliESDtrackCuts::GetTPCOnlyTrack(fEsdEv,etrack->GetID());
 110       if (!ntrack)
 111         continue;
 112       if (ntrack->Pt()<=0) {
 113         delete ntrack;
 114         continue;
 115       }
 116
 117       Double_t bfield[3] = {0,0,0};
 118       ntrack->GetBxByBz(bfield);
 119       AliExternalTrackParam exParam;
 120       Bool_t relate = ntrack->RelateToVertexBxByBz(vtxSPD,bfield,kVeryBig,&exParam);
 121       if (!relate) {
 122         delete ntrack;
 123         continue;
 124       }
 125
 126       // set the constraint parameters to the track
 127       ntrack->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
 128
 129       if (ntrack->Pt()<=0) {
 130         delete ntrack;
 131         continue;
 132       }
 133
 134       new ((*fTracks)[ntrnew++]) AliESDtrack(*ntrack);
 135       delete ntrack;
 136     }
 137
 138   } else {
 139
 140     am->LoadBranch("Tracks");
 141     Int_t ntr = fEsdEv->GetNumberOfTracks();
 142     for (Int_t i=0, ntrnew=0; i<ntr; ++i) {
 143       AliESDtrack *etrack = fEsdEv->GetTrack(i);
 144       if (!etrack)
 145         continue;
 146       if (fEsdTrackCuts->AcceptTrack(etrack)) {
 147         new ((*fTracks)[ntrnew]) AliESDtrack(*etrack);
 148         AliESDtrack *newTrack = static_cast<AliESDtrack*>(fTracks->At(ntrnew));
 149         newTrack->SetLabel(0);
 150         ++ntrnew;
 151       } else if (fHybridTrackCuts->AcceptTrack(etrack)) {
 152         if (etrack->GetConstrainedParam()) {
 153           new ((*fTracks)[ntrnew]) AliESDtrack(*etrack);
 154           AliESDtrack *newTrack = static_cast<AliESDtrack*>(fTracks->At(ntrnew));
 155           ++ntrnew;
 156           const AliExternalTrackParam* constrainParam = etrack->GetConstrainedParam();
 157           newTrack->Set(constrainParam->GetX(),
 158                         constrainParam->GetAlpha(),
 159                         constrainParam->GetParameter(),
 160                         constrainParam->GetCovariance());
 161           UInt_t status = etrack->GetStatus();
 162           if ((status&AliESDtrack::kITSrefit)==0)
 163             newTrack->SetLabel(2);
 164           else
 165             newTrack->SetLabel(1);
 166         }
 167       }
 168     }
 169   }
 170 }