PWG/EMCAL/AliEmcalAodTrackFilterTask.cxx

   1 // $Id: AliEmcalAodTrackFilterTask.cxx | Fri Dec 6 10:29:20 2013 +0100 | Constantin Loizides  $
   2 //
   3 // Class to filter Aod tracks
   4 //
   5 // Author: C.Loizides
   6
   7 #include <TClonesArray.h>
   8 #include <TRandom3.h>
   9 #include <AliAODEvent.h>
  10 #include <AliAODTrack.h>
  11 #include <AliAnalysisManager.h>
  12 #include <AliLog.h>
  13 #include <AliEMCALRecoUtils.h>
  14 #include "AliEmcalAodTrackFilterTask.h"
  15
  16 ClassImp(AliEmcalAodTrackFilterTask)
  17
  18 //________________________________________________________________________
  19 AliEmcalAodTrackFilterTask::AliEmcalAodTrackFilterTask() :
  20   AliAnalysisTaskSE("AliEmcalAodTrackFilterTask"),
  21   fTracksOutName("PicoTracks"),
  22   fTracksInName("tracks"),
  23   fMinTrackPt(0),
  24   fMaxTrackPt(1000),
  25   fMinTrackEta(-10),
  26   fMaxTrackEta(10),
  27   fMinTrackPhi(-10),
  28   fMaxTrackPhi(10),
  29   fTrackEfficiency(1),
  30   fIncludeNoITS(kTRUE),
  31   fUseNegativeLabels(kTRUE),
  32   fIsMC(kFALSE),
  33   fCutMaxFrShTPCClus(0.4),
  34   fModifyTrack(kTRUE),
  35   fDoPropagation(kFALSE),
  36   fDist(440),
  37   fTracksIn(0),
  38   fTracksOut(0)
  39 {
  40   // Constructor.
  41
  42   fAODfilterBits[0] = -1;
  43   fAODfilterBits[1] = -1;
  44 }
  45
  46 //________________________________________________________________________
  47 AliEmcalAodTrackFilterTask::AliEmcalAodTrackFilterTask(const char *name) :
  48   AliAnalysisTaskSE(name),
  49   fTracksOutName("PicoTracks"),
  50   fTracksInName("tracks"),
  51   fMinTrackPt(0),
  52   fMaxTrackPt(1000),
  53   fMinTrackEta(-10),
  54   fMaxTrackEta(10),
  55   fMinTrackPhi(-10),
  56   fMaxTrackPhi(10),
  57   fTrackEfficiency(1),
  58   fIncludeNoITS(kTRUE),
  59   fUseNegativeLabels(kTRUE),
  60   fIsMC(kFALSE),
  61   fCutMaxFrShTPCClus(0.4),
  62   fModifyTrack(kTRUE),
  63   fDoPropagation(kFALSE),
  64   fDist(440),
  65   fTracksIn(0),
  66   fTracksOut(0)
  67 {
  68   // Constructor.
  69
  70   fAODfilterBits[0] = -1;
  71   fAODfilterBits[1] = -1;
  72   fBranchNames = "AOD:tracks";
  73 }
  74
  75 //________________________________________________________________________
  76 AliEmcalAodTrackFilterTask::~AliEmcalAodTrackFilterTask()
  77 {
  78   // Destructor.
  79 }
  80
  81 //________________________________________________________________________
  82 void AliEmcalAodTrackFilterTask::UserCreateOutputObjects()
  83 {
  84   // Create my user objects.
  85
  86   fTracksOut = new TClonesArray("AliAODTrack");
  87   fTracksOut->SetName(fTracksOutName);
  88 }
  89
  90 //________________________________________________________________________
  91 void AliEmcalAodTrackFilterTask::UserExec(Option_t *)
  92 {
  93   // Main loop, called for each event.
  94
  95   AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
  96   if (!am) {
  97     AliError("Manager zero, returning");
  98     return;
  99   }
 100
 101   // retrieve tracks from input.
 102   if (!fTracksIn) {
 103     fTracksIn = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksInName));
 104     if (!fTracksIn) {
 105       AliError(Form("Could not retrieve tracks %s!", fTracksInName.Data()));
 106       return;
 107     }
 108     if (!fTracksIn->GetClass()->GetBaseClass("AliVParticle")) {
 109       AliError(Form("%s: Collection %s does not contain AliVParticle objects!", GetName(), fTracksInName.Data()));
 110       return;
 111     }
 112   }
 113
 114   // add tracks to event if not yet there
 115   fTracksOut->Delete();
 116   if (!(InputEvent()->FindListObject(fTracksOutName))) {
 117     InputEvent()->AddObject(fTracksOut);
 118   }
 119
 120   // loop over tracks
 121   const Int_t Ntracks = fTracksIn->GetEntriesFast();
 122   for (Int_t iTracks = 0, nacc = 0; iTracks < Ntracks; ++iTracks) {
 123
 124     AliAODTrack *track = static_cast<AliAODTrack*>(fTracksIn->At(iTracks));
 125
 126     if (!track)
 127       continue;
 128
 129     if (track->Pt() > fMaxTrackPt || track->Pt() < fMinTrackPt)
 130       continue;
 131
 132     if (track->Eta() < fMinTrackEta || track->Eta() > fMaxTrackEta ||
 133         track->Phi() < fMinTrackPhi || track->Phi() > fMaxTrackPhi)
 134       continue;
 135
 136     Int_t type = -1;
 137
 138     if (fAODfilterBits[0] < 0) {
 139       if (track->IsHybridGlobalConstrainedGlobal())
 140         type = 3;
 141       else /*not a good track*/
 142         continue;
 143     } else {
 144       if (track->TestFilterBit(fAODfilterBits[0])) {
 145         type = 0;
 146       } else if (track->TestFilterBit(fAODfilterBits[1])) {
 147         if ((track->GetStatus()&AliVTrack::kITSrefit)==0) {
 148           if (fIncludeNoITS)
 149             type = 2;
 150           else
 151             continue;
 152         } else {
 153           type = 1;
 154         }
 155       }
 156       else {/*not a good track*/
 157         continue;
 158       }
 159     }
 160
 161     if (fCutMaxFrShTPCClus > 0) {
 162       Double_t frac = Double_t(track->GetTPCnclsS()) / Double_t(track->GetTPCncls());
 163       if (frac > fCutMaxFrShTPCClus)
 164         continue;
 165     }
 166
 167     if (fTrackEfficiency < 1) {
 168       Double_t r = gRandom->Rndm();
 169       if (fTrackEfficiency < r)
 170         continue;
 171     }
 172
 173     Int_t label = 0;
 174     if (fIsMC) {
 175       if (fUseNegativeLabels)
 176         label = track->GetLabel();
 177       else
 178         label = TMath::Abs(track->GetLabel());
 179
 180       if (label == 0)
 181         AliDebug(2,Form("Track %d with label==0", iTracks));
 182     }
 183
 184     AliAODTrack *newt = new ((*fTracksOut)[nacc]) AliAODTrack(*track);
 185     if (fDoPropagation)
 186       AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface(newt,fDist);
 187     if (fModifyTrack) {
 188       newt->SetLabel(label);
 189       if (type==0) {
 190         newt->SetBit(BIT(27),0);
 191         newt->SetBit(BIT(28),0);
 192       } else if (type==1) {
 193         newt->SetBit(BIT(27),1);
 194         newt->SetBit(BIT(28),0);
 195       } else if (type==2) {
 196         newt->SetBit(BIT(27),0);
 197         newt->SetBit(BIT(28),1);
 198       } else if (type==3) {
 199         newt->SetBit(BIT(27),1);
 200         newt->SetBit(BIT(28),1);
 201       }
 202     }
 203     ++nacc;
 204   }
 205 }