PWG/EMCAL/AliAnalysisTaskEmcal.cxx

   1 // $Id: AliAnalysisTaskEmcal.cxx 56756 2012-05-30 05:03:02Z loizides $
   2 //
   3 // Emcal base analysis task.
   4 //
   5 // Author: S.Aiola
   6
   7 #include "AliAnalysisTaskEmcal.h"
   8
   9 #include <TChain.h>
  10 #include <TClonesArray.h>
  11 #include <TList.h>
  12 #include <TObject.h>
  13
  14 #include "AliAODEvent.h"
  15 #include "AliAnalysisManager.h"
  16 #include "AliCentrality.h"
  17 #include "AliEMCALGeometry.h"
  18 #include "AliESDEvent.h"
  19 #include "AliEmcalJet.h"
  20 #include "AliEmcalParticle.h"
  21 #include "AliEventplane.h"
  22 #include "AliInputEventHandler.h"
  23 #include "AliLog.h"
  24 #include "AliMCParticle.h"
  25 #include "AliVCluster.h"
  26 #include "AliVEventHandler.h"
  27 #include "AliVParticle.h"
  28
  29 ClassImp(AliAnalysisTaskEmcal)
  30
  31 //________________________________________________________________________
  32 AliAnalysisTaskEmcal::AliAnalysisTaskEmcal() :
  33   AliAnalysisTaskSE("AliAnalysisTaskEmcal"),
  34   fAnaType(kTPC),
  35   fInitialized(kFALSE),
  36   fCreateHisto(kTRUE),
  37   fTracksName(),
  38   fCaloName(),
  39   fMinCent(-999),
  40   fMaxCent(-999),
  41   fMinVz(-999),
  42   fMaxVz(-999),
  43   fOffTrigger(AliVEvent::kAny),
  44   fTrigClass(),
  45   fNbins(500),
  46   fMinBinPt(0),
  47   fMaxBinPt(250),
  48   fClusPtCut(0.15),
  49   fTrackPtCut(0.15),
  50   fMinTrackEta(-0.9),
  51   fMaxTrackEta(0.9),
  52   fMinTrackPhi(-10),
  53   fMaxTrackPhi(10),
  54   fClusTimeCutLow(-10),
  55   fClusTimeCutUp(10),
  56   fTracks(0),
  57   fCaloClusters(0),
  58   fCent(0),
  59   fCentBin(-1),
  60   fEPV0(-1.0),
  61   fEPV0A(-1.0),
  62   fEPV0C(-1.0),
  63   fNVertCont(0),
  64   fBeamType(kNA),
  65   fOutput(0)
  66 {
  67   // Default constructor.
  68
  69   fVertex[0] = 0;
  70   fVertex[1] = 0;
  71   fVertex[2] = 0;
  72 }
  73
  74 //________________________________________________________________________
  75 AliAnalysisTaskEmcal::AliAnalysisTaskEmcal(const char *name, Bool_t histo) :
  76   AliAnalysisTaskSE(name),
  77   fAnaType(kTPC),
  78   fInitialized(kFALSE),
  79   fCreateHisto(histo),
  80   fTracksName(),
  81   fCaloName(),
  82   fMinCent(-999),
  83   fMaxCent(-999),
  84   fMinVz(-999),
  85   fMaxVz(-999),
  86   fOffTrigger(AliVEvent::kAny),
  87   fTrigClass(),
  88   fNbins(500),
  89   fMinBinPt(0),
  90   fMaxBinPt(250),
  91   fClusPtCut(0.15),
  92   fTrackPtCut(0.15),
  93   fMinTrackEta(-0.9),
  94   fMaxTrackEta(0.9),
  95   fMinTrackPhi(-10),
  96   fMaxTrackPhi(10),
  97   fClusTimeCutLow(-10),
  98   fClusTimeCutUp(10),
  99   fTracks(0),
 100   fCaloClusters(0),
 101   fCent(0),
 102   fCentBin(-1),
 103   fEPV0(-1.0),
 104   fEPV0A(-1.0),
 105   fEPV0C(-1.0),
 106   fNVertCont(0),
 107   fBeamType(kNA),
 108   fOutput(0)
 109 {
 110   // Standard constructor.
 111
 112   fVertex[0] = 0;
 113   fVertex[1] = 0;
 114   fVertex[2] = 0;
 115
 116   if (fCreateHisto) {
 117     DefineOutput(1, TList::Class());
 118   }
 119 }
 120
 121 //________________________________________________________________________
 122 AliAnalysisTaskEmcal::~AliAnalysisTaskEmcal()
 123 {
 124   // Destructor
 125 }
 126
 127 //________________________________________________________________________
 128 void AliAnalysisTaskEmcal::UserExec(Option_t *)
 129 {
 130   // Main loop, called for each event.
 131
 132   if (!fInitialized)
 133     ExecOnce();
 134
 135   if (!fInitialized)
 136     return;
 137
 138   if (!RetrieveEventObjects())
 139     return;
 140
 141   if (!IsEventSelected())
 142     return;
 143
 144   if (!Run())
 145     return;
 146
 147   if (!FillHistograms())
 148     return;
 149
 150   if (fCreateHisto && fOutput) {
 151     // information for this iteration of the UserExec in the container
 152     PostData(1, fOutput);
 153   }
 154 }
 155
 156 //________________________________________________________________________
 157 Bool_t AliAnalysisTaskEmcal::AcceptCluster(AliVCluster *clus, Bool_t acceptMC) const
 158 {
 159   // Return true if cluster is accepted.
 160
 161   if (!clus)
 162     return kFALSE;
 163
 164   if (!clus->IsEMCAL())
 165     return kFALSE;
 166
 167   if (!acceptMC && clus->Chi2() == 100)
 168     return kFALSE;
 169
 170   if (clus->GetTOF() > fClusTimeCutUp || clus->GetTOF() < fClusTimeCutLow)
 171     return kFALSE;
 172
 173   TLorentzVector nPart;
 174   clus->GetMomentum(nPart, const_cast<Double_t*>(fVertex));
 175
 176   if (nPart.Et() < fClusPtCut)
 177     return kFALSE;
 178
 179   return kTRUE;
 180 }
 181
 182 //________________________________________________________________________
 183 Bool_t AliAnalysisTaskEmcal::AcceptEmcalPart(AliEmcalParticle *part, Bool_t acceptMC) const
 184 {
 185   // Return true if EMCal particle is accepted.
 186
 187   if (!part)
 188     return kFALSE;
 189
 190   if (fAnaType == kEMCAL && !part->IsEMCAL())
 191     return kFALSE;
 192
 193   if ((part->IsTrack() && part->Pt() < fTrackPtCut) || (part->IsCluster() && part->Pt() < fClusPtCut))
 194     return kFALSE;
 195
 196   if (!acceptMC && part->IsMC())
 197     return kFALSE;
 198
 199   return kTRUE;
 200 }
 201
 202 //________________________________________________________________________
 203 Bool_t AliAnalysisTaskEmcal::AcceptTrack(AliVTrack *track, Bool_t acceptMC) const
 204 {
 205   // Return true if track is accepted.
 206
 207   if (!track)
 208     return kFALSE;
 209
 210   if (!acceptMC && track->GetLabel() == 100)
 211     return kFALSE;
 212
 213   if (track->Pt() < fTrackPtCut)
 214     return kFALSE;
 215
 216   if (track->Eta() < fMinTrackEta || track->Eta() > fMaxTrackEta ||
 217       track->Phi() < fMinTrackPhi || track->Phi() > fMaxTrackPhi)
 218     return kFALSE;
 219
 220   return kTRUE;
 221 }
 222
 223 //________________________________________________________________________
 224 void AliAnalysisTaskEmcal::ExecOnce()
 225 {
 226   // Init the analysis.
 227
 228   if (!InputEvent()) {
 229     AliError(Form("%s: Could not retrieve event! Returning!", GetName()));
 230     return;
 231   }
 232
 233   if (!fCaloName.IsNull() && (fAnaType == kEMCAL || fAnaType == kEMCALOnly) && !fCaloClusters) {
 234     fCaloClusters =  dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCaloName));
 235     if (!fCaloClusters) {
 236       AliError(Form("%s: Could not retrieve clusters %s!", GetName(), fCaloName.Data()));
 237       return;
 238     } else {
 239       TClass *cl = fCaloClusters->GetClass();
 240       if (!cl->GetBaseClass("AliVCluster") && !cl->GetBaseClass("AliEmcalParticle")) {
 241         AliError(Form("%s: Collection %s does not contain AliVCluster nor AliEmcalParticle objects!", GetName(), fCaloName.Data()));
 242         fCaloClusters = 0;
 243         return;
 244       }
 245     }
 246   }
 247
 248   if (!fTracksName.IsNull() && fAnaType != kEMCALOnly && !fTracks) {
 249     fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName));
 250     if (!fTracks) {
 251       AliError(Form("%s: Could not retrieve tracks %s!", GetName(), fTracksName.Data()));
 252       return;
 253     } else {
 254       TClass *cl = fTracks->GetClass();
 255       if (!cl->GetBaseClass("AliVParticle") && !cl->GetBaseClass("AliEmcalParticle")) {
 256         AliError(Form("%s: Collection %s does not contain AliVParticle nor AliEmcalParticle objects!", GetName(), fTracksName.Data()));
 257         fTracks = 0;
 258         return;
 259       }
 260     }
 261   }
 262   SetInitialized();
 263 }
 264
 265 //_____________________________________________________
 266 AliAnalysisTaskEmcal::BeamType AliAnalysisTaskEmcal::GetBeamType()
 267 {
 268   // Get beam type : pp-AA-pA
 269   // ESDs have it directly, AODs get it from hardcoded run number ranges
 270
 271   AliESDEvent *esd = dynamic_cast<AliESDEvent*>(InputEvent());
 272   if (esd) {
 273     const AliESDRun *run = esd->GetESDRun();
 274     TString beamType = run->GetBeamType();
 275     if (beamType == "p-p")
 276       return kpp;
 277     else if (beamType == "A-A")
 278       return kAA;
 279     else if (beamType == "p-A")
 280       return kpA;
 281     else
 282       return kNA;
 283   } else {
 284     Int_t runNumber = InputEvent()->GetRunNumber();
 285     if ((runNumber >= 136851 && runNumber <= 139517) ||  // LHC10h
 286         (runNumber >= 166529 && runNumber <= 170593))    // LHC11h
 287     {
 288       return kAA;
 289     } else {
 290       return kpp;
 291     }
 292   }
 293 }
 294
 295 //________________________________________________________________________
 296 Bool_t AliAnalysisTaskEmcal::IsEventSelected()
 297 {
 298   // Check if event is selected
 299
 300   if (fOffTrigger != AliVEvent::kAny) {
 301     UInt_t res = 0;
 302     const AliESDEvent *eev = dynamic_cast<const AliESDEvent*>(InputEvent());
 303     if (eev) {
 304       res = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
 305     } else {
 306       const AliAODEvent *aev = dynamic_cast<const AliAODEvent*>(InputEvent());
 307       if (aev) {
 308         res = aev->GetHeader()->GetOfflineTrigger();
 309       }
 310     }
 311     if ((res & fOffTrigger) == 0)
 312       return kFALSE;
 313   }
 314
 315   if (!fTrigClass.IsNull()) {
 316     TString fired;
 317     const AliESDEvent *eev = dynamic_cast<const AliESDEvent*>(InputEvent());
 318     if (eev) {
 319       fired = eev->GetFiredTriggerClasses();
 320     } else {
 321       const AliAODEvent *aev = dynamic_cast<const AliAODEvent*>(InputEvent());
 322       if (aev) {
 323         fired = aev->GetFiredTriggerClasses();
 324       }
 325     }
 326     if (!fired.Contains("-B-"))
 327       return kFALSE;
 328     TObjArray *arr = fTrigClass.Tokenize("|");
 329     if (!arr)
 330       return kFALSE;
 331     Bool_t match = 0;
 332     for (Int_t i=0;i<arr->GetEntriesFast();++i) {
 333       TObject *obj = arr->At(i);
 334       if (!obj)
 335         continue;
 336       if (fired.Contains(obj->GetName())) {
 337         match = 1;
 338         break;
 339       }
 340     }
 341     delete arr;
 342     if (!match)
 343       return kFALSE;
 344   }
 345
 346   if ((fMinCent != -999) && (fMaxCent != -999)) {
 347     if (fCent<fMinCent)
 348       return kFALSE;
 349     if (fCent>fMaxCent)
 350       return kFALSE;
 351   }
 352
 353   if ((fMinVz != -999) && (fMaxVz != -999)) {
 354     if (fNVertCont == 0 )
 355       return kFALSE;
 356     Double_t vz = fVertex[2];
 357     if (vz<fMinVz)
 358       return kFALSE;
 359     if (vz>fMaxVz)
 360       return kFALSE;
 361   }
 362
 363   return kTRUE;
 364 }
 365
 366 //________________________________________________________________________
 367 TClonesArray *AliAnalysisTaskEmcal::GetArrayFromEvent(const char *name, const char *clname)
 368 {
 369   // Get array from event.
 370
 371   TClonesArray *arr = 0;
 372   TString sname(name);
 373   if (!sname.IsNull()) {
 374     arr = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(sname));
 375     if (!arr) {
 376       AliWarning(Form("%s: Could not retrieve array with name %s!", GetName(), name));
 377       return 0;
 378     }
 379   } else {
 380     return 0;
 381   }
 382
 383   if (!clname)
 384     return arr;
 385
 386   TString objname(arr->GetClass()->GetName());
 387   TClass cls(objname);
 388   if (!cls.InheritsFrom(clname)) {
 389     AliWarning(Form("%s: Objects of type %s in %s are not inherited from %s!",
 390                     GetName(), cls.GetName(), name, clname));
 391     return 0;
 392   }
 393   return arr;
 394 }
 395
 396 //________________________________________________________________________
 397 Bool_t AliAnalysisTaskEmcal::RetrieveEventObjects()
 398 {
 399   // Retrieve objects from event.
 400
 401   fVertex[0] = 0;
 402   fVertex[1] = 0;
 403   fVertex[2] = 0;
 404   fNVertCont = 0;
 405
 406   const AliVVertex *vert = InputEvent()->GetPrimaryVertex();
 407   if (vert) {
 408     vert->GetXYZ(fVertex);
 409     fNVertCont = vert->GetNContributors();
 410   }
 411
 412   fBeamType = GetBeamType();
 413
 414   if (fBeamType == kAA) {
 415     AliCentrality *aliCent = InputEvent()->GetCentrality();
 416     if (aliCent) {
 417       fCent = aliCent->GetCentralityPercentile("V0M");
 418       if      (fCent >=  0 && fCent <   10) fCentBin = 0;
 419       else if (fCent >= 10 && fCent <   30) fCentBin = 1;
 420       else if (fCent >= 30 && fCent <   50) fCentBin = 2;
 421       else if (fCent >= 50 && fCent <= 100) fCentBin = 3;
 422       else {
 423         AliWarning(Form("%s: Negative centrality: %f. Assuming 99", GetName(), fCent));
 424         fCentBin = 3;
 425       }
 426     } else {
 427       AliWarning(Form("%s: Could not retrieve centrality information! Assuming 99", GetName()));
 428       fCentBin = 3;
 429     }
 430     AliEventplane *aliEP = InputEvent()->GetEventplane();
 431     if (aliEP) {
 432       fEPV0  = aliEP->GetEventplane("V0" ,InputEvent());
 433       fEPV0A = aliEP->GetEventplane("V0A",InputEvent());
 434       fEPV0C = aliEP->GetEventplane("V0C",InputEvent());
 435     } else {
 436       AliWarning(Form("%s: Could not retrieve event plane information!", GetName()));
 437     }
 438   } else {
 439     fCent = 99;
 440     fCentBin = 0;
 441   }
 442
 443   return kTRUE;
 444 }