PWG3/hfe/AliAnalysisTaskHFE.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 /* $Id$ */
  17
  18 //
  19 // The analysis task:
  20 // Filling an AliCFContainer with the quantities pt, eta and phi
  21 // for tracks which survivied the particle cuts (MC resp. ESD tracks)
  22 // Track selection is done using the AliHFE package
  23 //
  24 // Author:
  25 //  Raphaelle Bailhache <R.Bailhache@gsi.de>
  26 //  Markus Fasel <M.Fasel@gsi.de>
  27 //  Matus Kalisky <matus.kalisky@cern.ch>
  28 //  MinJung Kweon <minjung@physi.uni-heidelberg.de>
  29 //
  30 #include <TAxis.h>
  31 #include <TCanvas.h>
  32 #include <TChain.h>
  33 #include <TDirectory.h>
  34 #include <TFile.h>
  35 #include <TH3D.h>
  36 #include <TIterator.h>
  37 #include <TList.h>
  38 #include <TLegend.h>
  39 #include <TMath.h>
  40 #include <TObjArray.h>
  41 #include <TParticle.h>
  42 #include <TProfile.h>
  43 #include <TString.h>
  44 #include <TF1.h>
  45 #include <TTree.h>
  46
  47 #include "AliAODInputHandler.h"
  48 #include "AliAODMCParticle.h"
  49 #include "AliAODTrack.h"
  50 #include "AliCFContainer.h"
  51 #include "AliCFManager.h"
  52 #include "AliESDEvent.h"
  53 #include "AliESDInputHandler.h"
  54 #include "AliESDpid.h"
  55 #include "AliESDtrack.h"
  56 #include "AliCentrality.h"
  57 #include "AliLog.h"
  58 #include "AliAnalysisManager.h"
  59 #include "AliMCEvent.h"
  60 #include "AliMCEventHandler.h"
  61 #include "AliMCParticle.h"
  62 #include "AliPID.h"
  63 #include "AliStack.h"
  64 #include "AliTriggerAnalysis.h"
  65 #include "AliVVertex.h"
  66
  67 #include "AliHFEcollection.h"
  68 #include "AliHFEcontainer.h"
  69 #include "AliHFEcuts.h"
  70 #include "AliHFEelecbackground.h"
  71 #include "AliHFEmcQA.h"
  72 #include "AliHFEpairs.h"
  73 #include "AliHFEpid.h"
  74 #include "AliHFEpidQAmanager.h"
  75 #include "AliHFEpostAnalysis.h"
  76 #include "AliHFEsecVtxs.h"
  77 #include "AliHFEsecVtx.h"
  78 #include "AliHFEsignalCuts.h"
  79 #include "AliHFEtaggedTrackAnalysis.h"
  80 #include "AliHFEtools.h"
  81 #include "AliHFEvarManager.h"
  82 #include "AliAnalysisTaskHFE.h"
  83
  84 ClassImp(AliAnalysisTaskHFE)
  85
  86 //____________________________________________________________
  87 AliAnalysisTaskHFE::AliAnalysisTaskHFE():
  88   AliAnalysisTaskSE("PID efficiency Analysis")
  89   , fQAlevel(0)
  90   , fPlugins(0)
  91   , fFillSignalOnly(kTRUE)
  92   , fRemovePileUp(kFALSE)
  93   , fIdentifiedAsPileUp(kFALSE)
  94   , fIdentifiedAsOutInz(kFALSE)
  95   , fPassTheEventCut(kFALSE)
  96   , fCentralityF(99.0)
  97   , fBackGroundFactorsFunction(NULL)
  98   , fContainer(NULL)
  99   , fVarManager(NULL)
 100   , fSignalCuts(NULL)
 101   , fCFM(NULL)
 102   , fTriggerAnalysis(NULL)
 103   , fPID(NULL)
 104   , fPIDqa(NULL)
 105   , fPIDpreselect(NULL)
 106   , fCuts(NULL)
 107   , fTaggedTrackCuts(NULL)
 108   , fCleanTaggedTrack(kFALSE)
 109   , fVariablesTRDTaggedTrack(kFALSE)
 110   , fCutspreselect(NULL)
 111   , fSecVtx(NULL)
 112   , fElecBackGround(NULL)
 113   , fMCQA(NULL)
 114   , fTaggedTrackAnalysis(NULL)
 115   , fQA(NULL)
 116   , fOutput(NULL)
 117   , fHistMCQA(NULL)
 118   , fHistSECVTX(NULL)
 119   , fHistELECBACKGROUND(NULL)
 120   , fQACollection(NULL)
 121 {
 122   //
 123   // Dummy constructor
 124   //
 125 }
 126
 127 //____________________________________________________________
 128 AliAnalysisTaskHFE::AliAnalysisTaskHFE(const char * name):
 129   AliAnalysisTaskSE(name)
 130   , fQAlevel(0)
 131   , fPlugins(0)
 132   , fFillSignalOnly(kTRUE)
 133   , fRemovePileUp(kFALSE)
 134   , fIdentifiedAsPileUp(kFALSE)
 135   , fIdentifiedAsOutInz(kFALSE)
 136   , fPassTheEventCut(kFALSE)
 137   , fCentralityF(99.0)
 138   , fBackGroundFactorsFunction(NULL)
 139   , fContainer(NULL)
 140   , fVarManager(NULL)
 141   , fSignalCuts(NULL)
 142   , fCFM(NULL)
 143   , fTriggerAnalysis(NULL)
 144   , fPID(NULL)
 145   , fPIDqa(NULL)
 146   , fPIDpreselect(NULL)
 147   , fCuts(NULL)
 148   , fTaggedTrackCuts(NULL)
 149   , fCleanTaggedTrack(kFALSE)
 150   , fVariablesTRDTaggedTrack(kFALSE)
 151   , fCutspreselect(NULL)
 152   , fSecVtx(NULL)
 153   , fElecBackGround(NULL)
 154   , fMCQA(NULL)
 155   , fTaggedTrackAnalysis(NULL)
 156   , fQA(NULL)
 157   , fOutput(NULL)
 158   , fHistMCQA(NULL)
 159   , fHistSECVTX(NULL)
 160   , fHistELECBACKGROUND(NULL)
 161   , fQACollection(0x0)
 162 {
 163   //
 164   // Default constructor
 165   //
 166   DefineOutput(1, TList::Class());
 167   DefineOutput(2, TList::Class());
 168
 169   fPID = new AliHFEpid("hfePid");
 170   fPIDqa = new AliHFEpidQAmanager;
 171   fVarManager = new AliHFEvarManager("hfeVarManager");
 172 }
 173
 174 //____________________________________________________________
 175 AliAnalysisTaskHFE::AliAnalysisTaskHFE(const AliAnalysisTaskHFE &ref):
 176   AliAnalysisTaskSE(ref)
 177   , fQAlevel(0)
 178   , fPlugins(0)
 179   , fFillSignalOnly(ref.fFillSignalOnly)
 180   , fRemovePileUp(ref.fRemovePileUp)
 181   , fIdentifiedAsPileUp(ref.fIdentifiedAsPileUp)
 182   , fIdentifiedAsOutInz(ref.fIdentifiedAsOutInz)
 183   , fPassTheEventCut(ref.fPassTheEventCut)
 184   , fCentralityF(ref.fCentralityF)
 185   , fBackGroundFactorsFunction(NULL)
 186   , fContainer(NULL)
 187   , fVarManager(NULL)
 188   , fSignalCuts(NULL)
 189   , fCFM(NULL)
 190   , fTriggerAnalysis(NULL)
 191   , fPID(NULL)
 192   , fPIDqa(NULL)
 193   , fPIDpreselect(NULL)
 194   , fCuts(NULL)
 195   , fTaggedTrackCuts(NULL)
 196   , fCleanTaggedTrack(ref.fCleanTaggedTrack)
 197   , fVariablesTRDTaggedTrack(ref.fVariablesTRDTaggedTrack)
 198   , fCutspreselect(NULL)
 199   , fSecVtx(NULL)
 200   , fElecBackGround(NULL)
 201   , fMCQA(NULL)
 202   , fTaggedTrackAnalysis(NULL)
 203   , fQA(NULL)
 204   , fOutput(NULL)
 205   , fHistMCQA(NULL)
 206   , fHistSECVTX(NULL)
 207   , fHistELECBACKGROUND(NULL)
 208   , fQACollection(NULL)
 209 {
 210   //
 211   // Copy Constructor
 212   //
 213   ref.Copy(*this);
 214 }
 215
 216 //____________________________________________________________
 217 AliAnalysisTaskHFE &AliAnalysisTaskHFE::operator=(const AliAnalysisTaskHFE &ref){
 218   //
 219   // Assignment operator
 220   //
 221   if(this == &ref)
 222     ref.Copy(*this);
 223   return *this;
 224 }
 225
 226 //____________________________________________________________
 227 void AliAnalysisTaskHFE::Copy(TObject &o) const {
 228   //
 229   // Copy into object o
 230   //
 231   AliAnalysisTaskHFE &target = dynamic_cast<AliAnalysisTaskHFE &>(o);
 232   target.fQAlevel = fQAlevel;
 233   target.fPlugins = fPlugins;
 234   target.fFillSignalOnly = fFillSignalOnly;
 235   target.fRemovePileUp = fRemovePileUp;
 236   target.fIdentifiedAsPileUp = fIdentifiedAsPileUp;
 237   target.fIdentifiedAsOutInz = fIdentifiedAsOutInz;
 238   target.fPassTheEventCut = fPassTheEventCut;
 239   target.fCentralityF = fCentralityF;
 240   target.fBackGroundFactorsFunction = fBackGroundFactorsFunction;
 241   target.fContainer = fContainer;
 242   target.fVarManager = fVarManager;
 243   target.fSignalCuts = fSignalCuts;
 244   target.fCFM = fCFM;
 245   target.fTriggerAnalysis = fTriggerAnalysis;
 246   target.fPID = fPID;
 247   target.fPIDqa = fPIDqa;
 248   target.fPIDpreselect = fPIDpreselect;
 249   target.fCuts = fCuts;
 250   target.fTaggedTrackCuts = fTaggedTrackCuts;
 251   target.fCleanTaggedTrack = fCleanTaggedTrack;
 252   target.fVariablesTRDTaggedTrack = fVariablesTRDTaggedTrack;
 253   target.fCutspreselect = fCutspreselect;
 254   target.fSecVtx = fSecVtx;
 255   target.fElecBackGround = fElecBackGround;
 256   target.fMCQA = fMCQA;
 257   target.fTaggedTrackAnalysis = fTaggedTrackAnalysis;
 258   target.fQA = fQA;
 259   target.fOutput = fOutput;
 260   target.fHistMCQA = fHistMCQA;
 261   target.fHistSECVTX = fHistSECVTX;
 262   target.fHistELECBACKGROUND = fHistELECBACKGROUND;
 263   target.fQACollection = fQACollection;
 264 }
 265
 266 //____________________________________________________________
 267 AliAnalysisTaskHFE::~AliAnalysisTaskHFE(){
 268   //
 269   // Destructor
 270   //
 271   if(fPID) delete fPID;
 272   if(fVarManager) delete fVarManager;
 273   if(fPIDqa) delete fPIDqa;
 274   if(fSignalCuts) delete fSignalCuts;
 275   if(fCFM) delete fCFM;
 276   if(fSecVtx) delete fSecVtx;
 277   if(fMCQA) delete fMCQA;
 278   if(fElecBackGround) delete fElecBackGround;
 279   if(fTriggerAnalysis) delete fTriggerAnalysis;
 280   if(fPIDpreselect) delete fPIDpreselect;
 281   if(fQA) delete fQA;
 282   if(fOutput) delete fOutput;
 283 }
 284
 285 //____________________________________________________________
 286 void AliAnalysisTaskHFE::UserCreateOutputObjects(){
 287   //
 288   // Creating output container and output objects
 289   // Here we also Initialize the correction framework container and
 290   // the objects for
 291   // - PID
 292   // - MC QA
 293   // - SecVtx
 294   // QA histograms are created if requested
 295   // Called once per worker
 296   //
 297   AliDebug(3, "Creating Output Objects");
 298   // Automatic determination of the analysis mode
 299   AliVEventHandler *inputHandler = dynamic_cast<AliVEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
 300   if(!TString(inputHandler->IsA()->GetName()).CompareTo("AliAODInputHandler")){
 301     SetAODAnalysis();
 302   } else {
 303     SetESDAnalysis();
 304     if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())
 305       SetHasMCData();
 306   }
 307   printf("Analysis Mode: %s Analysis\n", IsAODanalysis() ? "AOD" : "ESD");
 308   printf("MC Data available %s\n", HasMCData() ? "Yes" : "No");
 309
 310   // Enable Trigger Analysis
 311   fTriggerAnalysis = new AliTriggerAnalysis;
 312   fTriggerAnalysis->EnableHistograms();
 313   fTriggerAnalysis->SetAnalyzeMC(HasMCData());
 314
 315
 316   // Make lists for Output
 317   if(!fQA) fQA = new TList;
 318   fQA->SetOwner();
 319   if(!fOutput) fOutput = new TList;
 320   fOutput->SetOwner();
 321
 322   // First Part: Make QA histograms
 323   fQACollection = new AliHFEcollection("TaskQA", "QA histos from the Electron Task");
 324   fQACollection->CreateTH1F("nElectronTracksEvent", "Number of Electron Candidates", 100, 0, 100);
 325   fQACollection->CreateProfile("conr", "Electron PID contamination", 20, 0, 20);
 326   fQACollection->CreateTH1F("alpha_rec", "Alpha from reconstructed tracks with TRD hits", 36, -TMath::Pi(), TMath::Pi());
 327   fQACollection->CreateTH1F("alpha_sim", "Alpha from simulated electron tracks", 36, -TMath::Pi(), TMath::Pi());
 328   fQACollection->CreateTH1F("nElectron", "Number of electrons", 100, 0, 100);
 329   fQACollection->CreateProfile("pidquality", "TRD PID quality as function of momentum", 20, 0, 20);
 330   fQACollection->CreateProfile("ntrdclusters", "Number of TRD clusters as function of momentum", 20, 0, 20);
 331   fQACollection->CreateTH1F("chi2TRD","#chi2 per TRD cluster", 20, 0, 20);
 332   fQACollection->CreateTH1F("mccharge", "MC Charge", 200, -100, 100);
 333   fQACollection->CreateTH2F("radius", "Production Vertex", 100, 0.0, 5.0, 100, 0.0, 5.0);
 334   // Temporary histograms for TPC number of clusters for all signal tracks (MC true electrons) and for selected tracks (Markus Fasel)
 335   fQACollection->CreateTH2F("TPCclusters2_1_Signal", "TPCclusterInfo for findable clusters for 2 neighbors for signal tracks", 30, 0.1, 10., 162, 0., 161.);
 336   fQACollection->CreateTH2F("TPCclusters2_0_Signal", "TPCclusterInfo for the ratio for 2 neighbors for signal tracks", 30, 0.1, 10., 100, 0., 1.);
 337   fQACollection->CreateTH2F("TPCclusters2_1_Selected", "TPCclusterInfo for findable clusters for 2 neighbors for selected tracks", 30, 0.1, 10., 162, 0., 161.);
 338   fQACollection->CreateTH2F("TPCclusters2_0_Selected", "TPCclusterInfo for the ratio for 2 neighbors for selected tracks", 30, 0.1, 10., 110, 0., 1.1);
 339   fQACollection->CreateTH2F("TPCncls_Signal", "TPC Number of clusters for signal tracks", 30, 0.1, 10., 162, 0., 161.);
 340   fQACollection->CreateTH2F("TPCclr_Signal", "TPC cluster ratio for signal tracks", 30, 0.1, 10., 110, 0., 1.1);
 341   fQACollection->BinLogAxis("TPCclusters2_1_Signal", 0);
 342   fQACollection->BinLogAxis("TPCclusters2_0_Signal", 0);
 343   fQACollection->BinLogAxis("TPCclusters2_1_Selected", 0);
 344   fQACollection->BinLogAxis("TPCclusters2_0_Selected", 0);
 345   fQACollection->BinLogAxis("TPCncls_Signal", 0);
 346   fQACollection->BinLogAxis("TPCclr_Signal", 0);
 347
 348   InitPIDperformanceQA();
 349   InitContaminationQA();
 350   fQA->Add(fQACollection->GetList());
 351
 352   // Initialize PID
 353   fPID->SetHasMCData(HasMCData());
 354   if(!fPID->GetNumberOfPIDdetectors()) fPID->AddDetector("TPC", 0);
 355   fPID->InitializePID();
 356   if(IsQAOn(kPIDqa)){
 357     AliInfo("PID QA switched on");
 358     fPIDqa->Initialize(fPID);
 359     fQA->Add(fPIDqa->MakeList("HFEpidQA"));
 360   }
 361
 362   // Initialize correction Framework and Cuts
 363   const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
 364   fCFM = new AliCFManager;
 365   fCFM->SetNStepParticle(kNcutSteps);
 366   MakeParticleContainer();
 367   MakeEventContainer();
 368   // Temporary fix: Initialize particle cuts with NULL
 369   for(Int_t istep = 0; istep < kNcutSteps; istep++)
 370     fCFM->SetParticleCutsList(istep, NULL);
 371   if(!fCuts){
 372     AliWarning("Cuts not available. Default cuts will be used");
 373     fCuts = new AliHFEcuts;
 374     fCuts->CreateStandardCuts();
 375   }
 376   if(IsAODanalysis()) fCuts->SetAOD();
 377   // Make clone for V0 tagging step
 378   fCuts->Initialize(fCFM);
 379   if(fCuts->IsQAOn()) fQA->Add(fCuts->GetQAhistograms());
 380   fSignalCuts = new AliHFEsignalCuts("HFEsignalCuts", "HFE MC Signal definition");
 381   fVarManager->SetSignalCuts(fSignalCuts);
 382
 383   // add output objects to the List
 384   fOutput->AddAt(fContainer, 0);
 385   fOutput->AddAt(fCFM->GetEventContainer(), 1);
 386
 387   // mcQA----------------------------------
 388   if (HasMCData() && IsQAOn(kMCqa)) {
 389     AliInfo("MC QA on");
 390     if(!fMCQA) fMCQA = new AliHFEmcQA;
 391     if(!fHistMCQA) fHistMCQA = new TList();
 392     fHistMCQA->SetOwner();
 393     fMCQA->CreatDefaultHistograms(fHistMCQA);
 394     fQA->Add(fHistMCQA);
 395   }
 396
 397   // secvtx----------------------------------
 398   if (GetPlugin(kSecVtx)) {
 399     AliInfo("Secondary Vertex Analysis on");
 400     if(!fSecVtx) fSecVtx = new AliHFEsecVtx;
 401     fSecVtx->SetHasMCData(HasMCData());
 402
 403     if(!fHistSECVTX) fHistSECVTX = new TList();
 404     fHistSECVTX->SetOwner();
 405     fSecVtx->CreateHistograms(fHistSECVTX);
 406     fOutput->Add(fHistSECVTX);
 407   }
 408
 409   // background----------------------------------
 410   if (GetPlugin(kIsElecBackGround)) {
 411     AliInfo("Electron BackGround Analysis on");
 412     if(!fElecBackGround){
 413       AliWarning("ElecBackGround not available. Default elecbackground will be used");
 414       fElecBackGround = new AliHFEelecbackground;
 415     }
 416     fElecBackGround->SetHasMCData(HasMCData());
 417
 418     if(!fHistELECBACKGROUND) fHistELECBACKGROUND = new TList();
 419     fHistELECBACKGROUND->SetOwner();
 420     fElecBackGround->CreateHistograms(fHistELECBACKGROUND);
 421     fOutput->Add(fHistELECBACKGROUND);
 422   }
 423
 424   // tagged tracks
 425   if(GetPlugin(kTaggedTrackAnalysis)){
 426     AliInfo("Analysis on V0-tagged tracks enabled");
 427     fTaggedTrackAnalysis = new AliHFEtaggedTrackAnalysis;
 428     fTaggedTrackAnalysis->SetCuts(fTaggedTrackCuts);
 429     fTaggedTrackAnalysis->SetClean(fCleanTaggedTrack);
 430     if(fPIDqa->HasHighResolutionHistos())
 431       fTaggedTrackAnalysis->GetPIDqa()->SetHighResolutionHistos();
 432     fTaggedTrackAnalysis->SetPID(fPID);
 433     fTaggedTrackAnalysis->SetVariablesTRD(fVariablesTRDTaggedTrack);
 434     fTaggedTrackAnalysis->InitContainer();
 435     fOutput->Add(fTaggedTrackAnalysis->GetContainer());
 436     fQA->Add(fTaggedTrackAnalysis->GetPIDQA());
 437     fQA->Add(fTaggedTrackAnalysis->GetCutQA());
 438     fQA->Add(fTaggedTrackAnalysis->GetQAcollection());
 439   }
 440   PrintStatus();
 441 }
 442
 443 //____________________________________________________________
 444 void AliAnalysisTaskHFE::UserExec(Option_t *){
 445   //
 446   // Run the analysis
 447   //
 448   AliDebug(3, "Starting Single Event Analysis");
 449   if(!fInputEvent){
 450     AliError("Reconstructed Event not available");
 451     return;
 452   }
 453   if(HasMCData()){
 454     AliDebug(4, Form("MC Event: %p", fMCEvent));
 455     if(!fMCEvent){
 456       AliError("No MC Event, but MC Data required");
 457       return;
 458     }
 459   }
 460   if(!fCuts){
 461     AliError("HFE cuts not available");
 462     return;
 463   }
 464
 465   if(IsESDanalysis() && HasMCData()){
 466     // Protect against missing MC trees
 467     AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
 468     if(!mcH){
 469       AliError("No MC Event Handler available");
 470       return;
 471     }
 472     if(!mcH->InitOk()) return;
 473     if(!mcH->TreeK()) return;
 474     if(!mcH->TreeTR()) return;
 475   }
 476
 477   // need the centrality for everything (MC also)
 478   fCentralityF = 99.0;
 479   ReadCentrality();
 480   //printf("fCentralityF %f\n",fCentralityF);
 481
 482   // See if pile up and z in the range
 483   RejectionPileUpVertexRangeEventCut();
 484
 485   // Protect agains missing
 486   if(HasMCData()){
 487     //printf("Has MC data\n");
 488     fSignalCuts->SetMCEvent(fMCEvent);
 489     ProcessMC();  // Run the MC loop + MC QA in case MC Data are available
 490   }
 491
 492   if(IsAODanalysis()){
 493     AliAODpidUtil *aodworkingpid = AliHFEtools::GetDefaultAODPID(HasMCData());
 494     fPID->SetAODpid(aodworkingpid);
 495     ProcessAOD();
 496   } else {
 497     AliESDInputHandler *inH = dynamic_cast<AliESDInputHandler *>(fInputHandler);
 498     if(!inH){
 499       AliError("No ESD Input handler available");
 500       return;
 501     }
 502     AliESDpid *workingPID = inH->GetESDpid();
 503     if(!workingPID){
 504       AliDebug(1, "Using default ESD PID");
 505       workingPID = AliHFEtools::GetDefaultPID(HasMCData());
 506     } else {
 507       AliDebug(1, "Using ESD PID from the input handler");
 508     }
 509     fPID->SetESDpid(workingPID);
 510     if(fPIDpreselect) fPIDpreselect->SetESDpid(workingPID);
 511
 512     ProcessESD();
 513   }
 514   // Done!!!
 515   PostData(1, fOutput);
 516   PostData(2, fQA);
 517 }
 518
 519 //____________________________________________________________
 520 void AliAnalysisTaskHFE::Terminate(Option_t *){
 521   //
 522   // Terminate not implemented at the moment
 523   //
 524   if(GetPlugin(kPostProcess)){
 525     fOutput = dynamic_cast<TList *>(GetOutputData(1));
 526     fQA = dynamic_cast<TList *>(GetOutputData(2));
 527     if(!fOutput){
 528       AliError("Results not available");
 529       return;
 530     }
 531     if(!fQA){
 532       AliError("QA output not available");
 533       return;
 534     }
 535     fContainer = dynamic_cast<AliHFEcontainer *>(fOutput->FindObject("trackContainer"));
 536     if(!fContainer){
 537       AliError("Track container not found");
 538       return;
 539     }
 540     AliHFEpostAnalysis postanalysis;
 541     postanalysis.SetTaskResults(fContainer);
 542     TList *qalist = dynamic_cast<TList *>(fQA->FindObject("list_TaskQA"));
 543     if(!qalist){
 544       AliError("QA List not found");
 545       return;
 546     }
 547     postanalysis.SetTaskQA(qalist);
 548     printf("Running post analysis\n");
 549     //if(HasMCData())
 550     postanalysis.DrawMCSignal2Background();
 551     postanalysis.DrawEfficiency();
 552     postanalysis.DrawPIDperformance();
 553     postanalysis.DrawCutEfficiency();
 554
 555     if (GetPlugin(kIsElecBackGround)) {
 556       AliHFEelecbackground elecBackGround;
 557       TList *oe = 0x0;
 558       if(!(oe = (TList*)dynamic_cast<TList *>(fOutput->FindObject("HFEelecbackground")))){
 559         return;
 560       }
 561       elecBackGround.Load(oe);
 562       elecBackGround.Plot();
 563       elecBackGround.PostProcess();
 564     }
 565   }
 566 }
 567 //_______________________________________________________________
 568 Bool_t AliAnalysisTaskHFE::IsEventInBinZero() {
 569   //
 570   //
 571   //
 572
 573   //printf("test in IsEventInBinZero\n");
 574   if(!fInputEvent){
 575     AliError("Reconstructed Event not available");
 576     return kFALSE;
 577   }
 578
 579   // check vertex
 580   const AliVVertex *vertex = fInputEvent->GetPrimaryVertex();
 581   if(!vertex) return kTRUE;
 582   //if(vertex) return kTRUE;
 583
 584   // check tracks
 585   if(fInputEvent->GetNumberOfTracks()<=0) return kTRUE;
 586   //if(fInputEvent->GetNumberOfTracks()>0) return kTRUE;
 587
 588
 589   return kFALSE;
 590
 591 }
 592 //____________________________________________________________
 593 void AliAnalysisTaskHFE::ProcessMC(){
 594   //
 595   // Runs the MC Loop (filling the container for the MC Cut Steps with the observables pt, eta and phi)
 596   // In case MC QA is on also MC QA loop is done
 597   //
 598   AliDebug(3, "Processing MC Information");
 599   Double_t eventContainer [3];
 600   eventContainer[0] = fMCEvent->GetPrimaryVertex()->GetZ();
 601   eventContainer[2] = fCentralityF;
 602   //printf("z position is %f\n",eventContainer[0]);
 603   //if(fCFM->CheckEventCuts(AliHFEcuts::kEventStepGenerated, fMCEvent))
 604   fCFM->GetEventContainer()->Fill(eventContainer,AliHFEcuts::kEventStepGenerated);
 605   Int_t nElectrons = 0;
 606   if(IsESDanalysis()){
 607     if (HasMCData() && IsQAOn(kMCqa)) {
 608       AliDebug(2, "Running MC QA");
 609
 610       if(fMCEvent->Stack()){
 611               fMCQA->SetMCEvent(fMCEvent);
 612         fMCQA->SetGenEventHeader(fMCEvent->GenEventHeader());
 613         fMCQA->Init();
 614
 615         // loop over all tracks for decayed electrons
 616         for (Int_t igen = 0; igen < fMCEvent->GetNumberOfTracks(); igen++){
 617           TParticle* mcpart = fMCEvent->Stack()->Particle(igen);
 618           if(!mcpart) continue;
 619           fMCQA->GetQuarkKine(mcpart, igen, AliHFEmcQA::kCharm);
 620           fMCQA->GetQuarkKine(mcpart, igen, AliHFEmcQA::kBeauty);
 621           fMCQA->GetHadronKine(mcpart, AliHFEmcQA::kCharm);
 622           fMCQA->GetHadronKine(mcpart, AliHFEmcQA::kBeauty);
 623           fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kCharm,  AliHFEmcQA::kElectronPDG, 0); // no accept cut
 624           fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kBeauty, AliHFEmcQA::kElectronPDG, 0); // no accept cut
 625           fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kOthers, AliHFEmcQA::kElectronPDG, 0); // no accept cut
 626           if (TMath::Abs(mcpart->Eta()) < 0.9) {
 627             fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kCharm,  AliHFEmcQA::kElectronPDG, 1); // accept |eta|<0.9
 628             fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kBeauty, AliHFEmcQA::kElectronPDG, 1); // accept |eta|<0.9
 629             fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kOthers, AliHFEmcQA::kElectronPDG, 1); // accept |eta|<0.9
 630           }
 631           if (TMath::Abs(AliHFEtools::GetRapidity(mcpart)) < 0.5) {
 632             fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kCharm,  AliHFEmcQA::kElectronPDG, 2); // accept |y|<0.5
 633             fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kBeauty, AliHFEmcQA::kElectronPDG, 2); // accept |y|<0.5
 634             fMCQA->GetDecayedKine(mcpart, AliHFEmcQA::kOthers, AliHFEmcQA::kElectronPDG, 2); // accept |y|<0.5
 635           }
 636         }
 637         fMCQA->EndOfEventAna(AliHFEmcQA::kCharm);
 638         fMCQA->EndOfEventAna(AliHFEmcQA::kBeauty);
 639       }
 640
 641     } // end of MC QA loop
 642     // -----------------------------------------------------------------
 643     fCFM->SetMCEventInfo(fMCEvent);
 644     // fCFM->CheckEventCuts(AliCFManager::kEvtRecCuts, fESD);
 645   } else {
 646     fCFM->SetMCEventInfo(fInputEvent);
 647   }
 648   // Run MC loop
 649   AliVParticle *mctrack = NULL;
 650   AliDebug(3, Form("Number of Tracks: %d", fMCEvent->GetNumberOfTracks()));
 651   for(Int_t imc = 0; imc <fMCEvent->GetNumberOfTracks(); imc++){
 652     if(!(mctrack = fMCEvent->GetTrack(imc))) continue;
 653     AliDebug(4, "Next MC Track");
 654     if(ProcessMCtrack(mctrack)) nElectrons++;
 655   }
 656
 657   // fCFM->CheckEventCuts(AliCFManager::kEvtRecCuts, fESD);
 658   fQACollection->Fill("nElectron", nElectrons);
 659 }
 660
 661 //____________________________________________________________
 662 void AliAnalysisTaskHFE::ProcessESD(){
 663   //
 664   // Run Analysis of reconstructed event in ESD Mode
 665   // Loop over Tracks, filter according cut steps defined in AliHFEcuts
 666   //
 667   AliDebug(3, "Processing ESD Event");
 668   AliESDEvent *fESD = dynamic_cast<AliESDEvent *>(fInputEvent);
 669   if(!fESD){
 670     AliError("ESD Event required for ESD Analysis")
 671     return;
 672   }
 673
 674   // Set magnetic field if V0 task on
 675   if(fTaggedTrackAnalysis) fTaggedTrackAnalysis->SetMagneticField(fESD->GetMagneticField());
 676
 677   // Do event Normalization
 678   Double_t eventContainer[3];
 679   eventContainer[0] = 0.0;
 680   if(fESD->GetPrimaryVertexTracks()) eventContainer[0] = fESD->GetPrimaryVertexTracks()->GetZ();
 681   eventContainer[1] = 0.;
 682   eventContainer[2] = fCentralityF;
 683   if(fTriggerAnalysis->IsOfflineTriggerFired(fESD, AliTriggerAnalysis::kV0AND))
 684     eventContainer[1] = 1.;
 685
 686   //
 687   fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepRecNoCut);
 688
 689   //
 690   if(fIdentifiedAsPileUp) return;
 691   fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepRecNoPileUp);
 692
 693   //
 694   if(fIdentifiedAsOutInz) return;
 695   fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepZRange);
 696
 697   //
 698   if(!fPassTheEventCut) return;
 699   fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepReconstructed);
 700
 701
 702
 703   fContainer->NewEvent();
 704
 705   if (GetPlugin(kIsElecBackGround)) {
 706     fElecBackGround->SetEvent(fESD);
 707   }
 708   if (GetPlugin(kSecVtx)) {
 709     fSecVtx->SetEvent(fESD);
 710     fSecVtx->GetPrimaryCondition();
 711   }
 712
 713   if(HasMCData()){
 714     if (GetPlugin(kSecVtx)) {
 715       fSecVtx->SetMCEvent(fMCEvent);
 716       fSecVtx->SetMCQA(fMCQA);
 717     }
 718     if (GetPlugin(kIsElecBackGround)) {
 719       fElecBackGround->SetMCEvent(fMCEvent);
 720     }
 721   }
 722
 723   Double_t container[10];
 724   memset(container, 0, sizeof(Double_t) * 10);
 725   // container for the output THnSparse
 726   Double_t dataE[6]; // [pT, eta, Phi, type, 'C' or 'B']
 727   Int_t nElectronCandidates = 0;
 728   AliESDtrack *track = NULL, *htrack = NULL;
 729   AliMCParticle *mctrack = NULL;
 730   TParticle* mctrack4QA = NULL;
 731   Int_t pid = 0;
 732
 733   Bool_t signal = kTRUE;
 734
 735   fCFM->SetRecEventInfo(fESD);
 736   // Electron background analysis
 737   if (GetPlugin(kIsElecBackGround)) {
 738
 739     AliDebug(2, "Running BackGround Analysis");
 740
 741     fElecBackGround->Reset();
 742
 743   } // end of electron background analysis
 744   //
 745   // Loop ESD
 746   //
 747   AliDebug(3, Form("Number of Tracks: %d", fESD->GetNumberOfTracks()));
 748   for(Int_t itrack = 0; itrack < fESD->GetNumberOfTracks(); itrack++){
 749     AliDebug(4, "New ESD track");
 750     track = fESD->GetTrack(itrack);
 751
 752     // fill counts of v0-identified particles
 753     Int_t v0pid = -1;
 754     if(track->TestBit(BIT(14))) v0pid = AliPID::kElectron;
 755     else if(track->TestBit(BIT(15))) v0pid = AliPID::kPion;
 756     else if(track->TestBit(BIT(16))) v0pid = AliPID::kProton;
 757     // here the tagged track analysis will run
 758     if(fTaggedTrackAnalysis && v0pid > -1){
 759       AliDebug(1, Form("Track identified as %s", AliPID::ParticleName(v0pid)));
 760       fTaggedTrackAnalysis->ProcessTrack(track, v0pid);
 761     }
 762
 763     AliDebug(3, Form("Doing track %d, %p", itrack, track));
 764
 765     //////////////////////////////////////
 766     // preselect
 767     /////////////////////////////////////
 768     if(fPIDpreselect && fCutspreselect) {
 769       if(!PreSelectTrack(track)) continue;
 770     }
 771
 772     signal = kTRUE;
 773
 774     // Fill step without any cut
 775
 776     if(HasMCData()){
 777       // Check if it is electrons near the vertex
 778       if(!(mctrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(track->GetLabel()))))) continue;
 779       mctrack4QA = mctrack->Particle();
 780
 781       if(fFillSignalOnly && !fCFM->CheckParticleCuts(AliHFEcuts::kStepMCGenerated, mctrack)) signal = kFALSE;
 782       else AliDebug(3, "Signal Electron");
 783     }
 784     // Cache new Track information inside the var manager
 785     fVarManager->NewTrack(track, mctrack, fCentralityF, -1, signal);
 786
 787     if(signal) {
 788       fVarManager->FillContainer(fContainer, "recTrackContReco", AliHFEcuts::kStepRecNoCut, kFALSE);
 789       fVarManager->FillContainer(fContainer, "recTrackContMC", AliHFEcuts::kStepRecNoCut, kTRUE);
 790       if((track->GetStatus() & AliESDtrack::kTPCout)
 791           && (TMath::Abs(track->Eta()) < 0.8)
 792           && (track->GetKinkIndex(0) == 0)){
 793         fQACollection->Fill("TPCclusters2_1_Signal", track->Pt(), track->GetTPCClusterInfo(2,1));
 794         fQACollection->Fill("TPCclusters2_0_Signal", track->Pt(), track->GetTPCNclsF() > 0 ?  track->GetTPCClusterInfo(2,1)/track->GetTPCNclsF() : 0.);
 795         fQACollection->Fill("TPCncls_Signal", track->Pt(), track->GetTPCNcls());
 796         fQACollection->Fill("TPCclr_Signal", track->Pt(), track->GetTPCNclsF() > 0 ? static_cast<Double_t>(track->GetTPCNcls())/static_cast<Double_t>(track->GetTPCNclsF()) : 0.);
 797       }
 798     }
 799
 800     // RecKine: ITSTPC cuts
 801     if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
 802
 803     // Check TRD criterions (outside the correction framework)
 804     if(track->GetTRDncls()){
 805       fQACollection->Fill("chi2TRD", track->GetTRDchi2()/track->GetTRDncls());
 806       fQACollection->Fill("alpha_rec", track->GetAlpha());
 807       fQACollection->Fill("pidquality", container[0], track->GetTRDpidQuality());
 808       fQACollection->Fill("ntrdclusters", container[0], track->GetTRDncls());
 809     }
 810
 811
 812     // RecPrim
 813     if(track->GetKinkIndex(0) != 0) continue; // Quick and dirty fix to reject both kink mothers and daughters
 814     if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
 815
 816     // HFEcuts: ITS layers cuts
 817     if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
 818
 819     if(HasMCData() && IsQAOn(kMCqa)) {
 820       // mc qa for after the reconstruction cuts
 821       AliDebug(2, "Running MC QA");
 822       fMCQA->GetDecayedKine(mctrack4QA, AliHFEmcQA::kCharm,  AliHFEmcQA::kElectronPDG, 3);  // charm
 823       fMCQA->GetDecayedKine(mctrack4QA, AliHFEmcQA::kBeauty,  AliHFEmcQA::kElectronPDG, 3); // beauty
 824       fMCQA->GetDecayedKine(mctrack4QA, AliHFEmcQA::kOthers,  AliHFEmcQA::kElectronPDG, 3); // beauty
 825     }
 826
 827     // HFEcuts: Nb of tracklets TRD0
 828     if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTRD, track)) continue;
 829
 830     // Fill correlation maps before PID
 831     if(signal && fContainer->GetCorrelationMatrix("correlationstepbeforePID")) {
 832       //printf("Fill correlation maps before PID\n");
 833       fVarManager->FillCorrelationMatrix(fContainer->GetCorrelationMatrix("correlationstepbeforePID"));
 834     }
 835
 836     if (HasMCData() && IsQAOn(kMCqa)) {
 837       // mc qa for after the reconstruction and pid cuts
 838       AliDebug(2, "Running MC QA");
 839       fMCQA->GetDecayedKine(mctrack4QA, AliHFEmcQA::kCharm,  AliHFEmcQA::kElectronPDG, 4);  // charm
 840       fMCQA->GetDecayedKine(mctrack4QA, AliHFEmcQA::kBeauty,  AliHFEmcQA::kElectronPDG, 4); // beauty
 841       fMCQA->GetDecayedKine(mctrack4QA, AliHFEmcQA::kOthers,  AliHFEmcQA::kElectronPDG, 4); // beauty
 842     }
 843
 844     if(HasMCData()){
 845       FillProductionVertex(track);
 846     }
 847
 848     // track accepted, do PID
 849     AliHFEpidObject hfetrack;
 850     hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
 851     hfetrack.SetRecTrack(track);
 852     if(HasMCData()) hfetrack.SetMCTrack(mctrack);
 853     hfetrack.SetCentrality(fCentralityF);
 854     fPID->SetVarManager(fVarManager);
 855     if(!fPID->IsSelected(&hfetrack, fContainer, "recTrackCont", fPIDqa)) continue;
 856     nElectronCandidates++;
 857     fQACollection->Fill("TPCclusters2_1_Selected", track->Pt(), track->GetTPCClusterInfo(2,1));
 858     fQACollection->Fill("TPCclusters2_0_Selected", track->Pt(), track->GetTPCClusterInfo(2,0));
 859
 860     // Fill Histogram for Hadronic Background
 861     if(HasMCData()){
 862       if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) != 11))
 863         fVarManager->FillContainer(fContainer, "hadronicBackground", UInt_t(0), kFALSE);
 864     }
 865
 866     // Fill Containers
 867     if(signal) {
 868       // Apply weight for background contamination
 869       if(fBackGroundFactorsFunction) {
 870               Double_t weightBackGround = fBackGroundFactorsFunction->Eval(TMath::Abs(track->P()));
 871               if(weightBackGround < 0.0) weightBackGround = 0.0;
 872         else if(weightBackGround > 1.0) weightBackGround = 1.0;
 873         // weightBackGround as special weight
 874         fVarManager->FillContainer(fContainer, "hadronicBackground", 1, kFALSE, weightBackGround);
 875       }
 876       fVarManager->FillCorrelationMatrix(fContainer->GetCorrelationMatrix("correlationstepafterPID"));
 877     }
 878
 879     if(GetPlugin(kSecVtx)) {
 880       AliDebug(2, "Running Secondary Vertex Analysis");
 881       fSecVtx->Process(track);
 882     }
 883
 884     if(HasMCData()){
 885       dataE[0] = track->Pt();
 886       dataE[1] = track->Eta();
 887       dataE[2] = track->Phi();
 888       dataE[3] = track->Charge();
 889       dataE[4] = -1.;
 890       dataE[5] = -1.;
 891
 892       // Track selected: distinguish between true and fake
 893       AliDebug(1, Form("Candidate Selected, filling THnSparse, PID: %d\n", mctrack->Particle()->GetPdgCode()));
 894       if((pid = TMath::Abs(mctrack->Particle()->GetPdgCode())) == 11){
 895         Int_t type = 0;
 896         if(fSignalCuts->IsCharmElectron(track))
 897           type = 1;
 898         else if(fSignalCuts->IsBeautyElectron(track))
 899           type = 2;
 900         AliDebug(1, Form("Type: %d\n", type));
 901         if(type){
 902                 dataE[5] = type; // beauty[1] or charm[2]
 903                 dataE[4] = 2;  // signal electron
 904         }
 905         else{
 906                 dataE[4] = 1; // not a signal electron
 907                 dataE[5] = 0;
 908         }
 909       }
 910       else {
 911         // Fill THnSparse with the information for Fake Electrons
 912         dataE[4] = 0;
 913         dataE[5] = 0;
 914       }
 915       // fill the performance THnSparse, if the mc origin could be defined
 916       if(dataE[4] > -1){
 917         AliDebug(1, Form("Entries: [%.3f|%.3f|%.3f|%f|%f|%f]\n", dataE[0],dataE[1],dataE[2],dataE[3],dataE[4],dataE[5]));
 918         fQACollection->Fill("PIDperformance", dataE);
 919       }
 920     }
 921     // Electron background analysis
 922     if (GetPlugin(kIsElecBackGround)) {
 923
 924       AliDebug(2, "Running BackGround Analysis");
 925
 926       for(Int_t jtrack = 0; jtrack < fESD->GetNumberOfTracks(); jtrack++){
 927               htrack = fESD->GetTrack(jtrack);
 928         if ( itrack == jtrack ) continue;
 929         fElecBackGround->PairAnalysis(track, htrack);
 930       }
 931     } // end of electron background analysis
 932
 933     if (GetPlugin(kDEstep)) {
 934       if(HasMCData()){
 935         if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) != 11)){
 936           fQACollection->Fill("hadronsBeforeIPcut",track->Pt());
 937           fQACollection->Fill("hadronsBeforeIPcutMC",mctrack->Pt());
 938         }
 939       }
 940       // Fill Containers for impact parameter analysis
 941       if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepHFEcutsDca + AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack,track)) continue;
 942       if(signal) {
 943         fVarManager->FillContainer(fContainer, "recTrackContDEReco", AliHFEcuts::kStepHFEcutsDca, kFALSE);
 944         fVarManager->FillContainer(fContainer, "recTrackContDEMC", AliHFEcuts::kStepHFEcutsDca, kTRUE);
 945         fVarManager->FillCorrelationMatrix(fContainer->GetCorrelationMatrix("correlationstepafterDE"));
 946       }
 947       if(HasMCData()){
 948         if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) != 11)){
 949           fVarManager->FillContainer(fContainer, "hadronicBackground", 2, kFALSE);
 950           fQACollection->Fill("hadronsAfterIPcut",track->Pt());
 951           fQACollection->Fill("hadronsAfterIPcutMC",mctrack->Pt());
 952         }
 953       }
 954     }
 955
 956   }
 957   fQACollection->Fill("nElectronTracksEvent", nElectronCandidates);
 958 }
 959
 960 //____________________________________________________________
 961 void AliAnalysisTaskHFE::ProcessAOD(){
 962   //
 963   // Run Analysis in AOD Mode
 964   // Function is still in development
 965   //
 966   AliDebug(3, "Processing AOD Event");
 967   Double_t eventContainer[3];
 968   eventContainer[0] = fInputEvent->GetPrimaryVertex()->GetZ();
 969   eventContainer[1] = 1.; // No Information available in AOD analysis, assume all events have V0AND
 970   eventContainer[2] = fCentralityF;
 971
 972   AliAODEvent *fAOD = dynamic_cast<AliAODEvent *>(fInputEvent);
 973   if(!fAOD){
 974     AliError("AOD Event required for AOD Analysis")
 975       return;
 976   }
 977
 978   //
 979   fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepRecNoCut);
 980
 981   //
 982   if(fIdentifiedAsPileUp) return;
 983   fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepRecNoPileUp);
 984
 985   //
 986   if(fIdentifiedAsOutInz) return;
 987   fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepZRange);
 988
 989   //
 990   if(!fPassTheEventCut) return;
 991   fCFM->GetEventContainer()->Fill(eventContainer, AliHFEcuts::kEventStepReconstructed);
 992
 993   fContainer->NewEvent();
 994
 995   AliAODTrack *track = NULL;
 996   AliAODMCParticle *mctrack = NULL;
 997   Double_t dataE[6]; // [pT, eta, Phi, Charge, type, 'C' or 'B']
 998   Int_t nElectronCandidates = 0;
 999   Int_t pid;
1000   Bool_t signal;
1001   for(Int_t itrack = 0; itrack < fAOD->GetNumberOfTracks(); itrack++){
1002     track = fAOD->GetTrack(itrack);
1003     if(!track) continue;
1004     if(track->GetFlags() != 1<<4) continue;  // Only process AOD tracks where the HFE is set
1005
1006     signal = kTRUE;
1007     if(HasMCData()){
1008
1009       Int_t label = TMath::Abs(track->GetLabel());
1010       if(label)
1011         mctrack = dynamic_cast<AliAODMCParticle *>(fMCEvent->GetTrack(label));
1012         if(fFillSignalOnly && !fCFM->CheckParticleCuts(AliHFEcuts::kStepMCGenerated, mctrack)) signal = kFALSE;
1013     }
1014     fVarManager->NewTrack(track, mctrack, fCentralityF, -1, kTRUE);
1015     // track accepted, do PID
1016     AliHFEpidObject hfetrack;
1017     hfetrack.SetAnalysisType(AliHFEpidObject::kAODanalysis);
1018     hfetrack.SetRecTrack(track);
1019     if(HasMCData()) hfetrack.SetMCTrack(mctrack);
1020     hfetrack.SetCentrality(fCentralityF);
1021     fPID->SetVarManager(fVarManager);
1022     if(!fPID->IsSelected(&hfetrack, fContainer, "recTrackCont", fPIDqa)) continue;    // we will do PID here as soon as possible
1023     // Apply weight for background contamination
1024     Double_t weightBackGround = 1.0;
1025     if(fBackGroundFactorsFunction) {
1026       weightBackGround = fBackGroundFactorsFunction->Eval(TMath::Abs(track->P()));
1027       if(weightBackGround < 0.0) weightBackGround = 0.0;
1028     }
1029     fVarManager->FillContainer(fContainer, "hadronicBackground", 1, kFALSE);
1030     nElectronCandidates++;
1031     if(HasMCData()){
1032       dataE[0] = track->Pt();
1033       dataE[1] = track->Eta();
1034       dataE[2] = track->Phi();
1035       dataE[3] = track->Charge();
1036       dataE[4] = -1;
1037       dataE[5] = -1;
1038      // Track selected: distinguish between true and fake
1039       AliDebug(1, Form("Candidate Selected, filling THnSparse, PID: %d\n", mctrack->GetPdgCode()));
1040       if((pid = TMath::Abs(mctrack->GetPdgCode())) == 11){
1041
1042         Int_t type = 0;
1043         if(fSignalCuts->IsCharmElectron(track))
1044           type = 1;
1045         else if(fSignalCuts->IsBeautyElectron(track))
1046           type = 2;
1047         AliDebug(1, Form("Type: %d\n", type));
1048         if(type){
1049                 dataE[5] = type; // beauty[1] or charm[2]
1050                 dataE[4] = 2;  // signal electron
1051         }
1052         else{
1053                 dataE[4] = 1; // not a signal electron
1054                 dataE[5] = 0;
1055         }
1056       }
1057       else {
1058         // Fill THnSparse with the information for Fake Electrons
1059         dataE[4] = 0;
1060         dataE[5] = 0;
1061       }
1062       // fill the performance THnSparse, if the mc origin could be defined
1063       if(dataE[4] > -1){
1064         AliDebug(1, Form("Entries: [%.3f|%.3f|%.3f|%f|%f|%f]\n", dataE[0],dataE[1],dataE[2],dataE[3],dataE[4],dataE[5]));
1065         fQACollection->Fill("PIDperformance", dataE);
1066       }
1067     }
1068   }
1069   fQACollection->Fill("nElectronTracksEvent", nElectronCandidates);
1070 }
1071
1072 //____________________________________________________________
1073 Bool_t AliAnalysisTaskHFE::ProcessMCtrack(AliVParticle *track){
1074   //
1075   // Filter the Monte Carlo Track
1076   // Additionally Fill a THnSparse for Signal To Background Studies
1077   // Works for AOD and MC analysis Type
1078   //
1079   fVarManager->NewTrack(track, NULL, -1, kTRUE);
1080   Double_t signalContainer[6];
1081
1082   signalContainer[0] = track->Pt();
1083   signalContainer[1] = track->Eta();
1084   signalContainer[2] = track->Phi();
1085   signalContainer[3] = track->Charge()/3;
1086
1087  Double_t vertex[3]; // Production vertex cut to mask gammas which are NOT supposed to have hits in the first ITS layer(s)
1088   if(IsESDanalysis()){
1089     AliMCParticle *mctrack = dynamic_cast<AliMCParticle *>(track);
1090     if(mctrack){
1091       vertex[0] = mctrack->Particle()->Vx();
1092       vertex[1] = mctrack->Particle()->Vy();
1093     }
1094   } else {
1095     AliAODMCParticle *aodmctrack = dynamic_cast<AliAODMCParticle *>(track);
1096     if(aodmctrack) aodmctrack->XvYvZv(vertex);
1097   }
1098
1099   if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepMCGenerated, track)) return kFALSE;
1100   fQACollection->Fill("mccharge", signalContainer[3]);
1101   fVarManager->FillContainer(fContainer, "MCTrackCont", AliHFEcuts::kStepMCGenerated, kFALSE);
1102   signalContainer[4] = 0;
1103   if(fSignalCuts->IsSelected(track)){
1104     //fVarManager->FillContainer(fContainer, "MCTrackCont", AliHFEcuts::kStepMCsignal, kFALSE);
1105     // Filling of the Signal/Background histogram using the
1106     // definition of the codes for charm and beauty as below in
1107     // th crearion of the histogram
1108     if(fSignalCuts->IsCharmElectron(track))
1109       signalContainer[4] = 1;
1110     else
1111       signalContainer[4] = 2;
1112   } else {
1113     signalContainer[4] = 0; // (and other background)
1114   }
1115   signalContainer[5] = 0;
1116   // apply cut on the sqrt of the production vertex
1117   Double_t radVertex = TMath::Sqrt(vertex[0]*vertex[0] + vertex[1] * vertex[1]);
1118   if(radVertex < 3.5){
1119     // Within first ITS layer(2) -> Background we cannot reject by ITS cut, let it pass
1120     signalContainer[5] = 1;
1121   } else if (radVertex < 7.5){
1122     signalContainer[5] = 2;
1123   }
1124   fQACollection->Fill("SignalToBackgroundMC", signalContainer);
1125   fQACollection->Fill("alpha_sim", track->Phi() - TMath::Pi());
1126
1127   // Step GeneratedZOutNoPileUp
1128   if((fIdentifiedAsPileUp) || (fIdentifiedAsOutInz)) return kFALSE;
1129   fVarManager->FillContainer(fContainer, "MCTrackCont", AliHFEcuts::kStepMCGeneratedZOutNoPileUp, kFALSE);
1130
1131   // Step Generated Event Cut
1132   if(!fPassTheEventCut) return kFALSE;
1133   fVarManager->FillContainer(fContainer, "MCTrackCont", AliHFEcuts::kStepMCGeneratedEventCut, kFALSE);
1134
1135   if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepMCInAcceptance, track)) return kFALSE;
1136   fVarManager->FillContainer(fContainer, "MCTrackCont", AliHFEcuts::kStepMCInAcceptance, kFALSE);
1137   return kTRUE;
1138 }
1139
1140 //____________________________________________________________
1141 Bool_t AliAnalysisTaskHFE::PreSelectTrack(AliESDtrack *track) const {
1142   //
1143   // Preselect tracks
1144   //
1145
1146
1147   Bool_t survived = kTRUE;
1148
1149   if(!fCutspreselect->CheckParticleCuts(AliHFEcuts::kStepRecKineITSTPC, track)) {
1150     survived = kFALSE;
1151     //printf("Did not pass AliHFEcuts::kStepRecKineITSTPC\n");
1152   }
1153   //else printf("Pass AliHFEcuts::kStepRecKineITSTPC\n");
1154   if(!fCutspreselect->CheckParticleCuts(AliHFEcuts::kStepRecPrim, track)) {
1155     survived = kFALSE;
1156     //printf("Did not pass AliHFEcuts::kStepRecPrim\n");
1157   }
1158   //else printf("Pass AliHFEcuts::kStepRecPrim\n");
1159   if(!fCutspreselect->CheckParticleCuts(AliHFEcuts::kStepHFEcutsITS, track)) {
1160     survived = kFALSE;
1161     //printf("Did not pass AliHFEcuts::kStepHFEcutsITS\n");
1162   }
1163   //else printf("Pass AliHFEcuts::kStepHFEcutsITS\n");
1164   if(!fCutspreselect->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTRD, track)) {
1165     survived = kFALSE;
1166     //printf("Did not pass AliHFEcuts::kStepHFEcutsTRD\n");
1167   }
1168   //else printf("Pass AliHFEcuts::kStepHFEcutsTRD\n");
1169
1170   if(survived){
1171     // Apply PID
1172     AliHFEpidObject hfetrack;
1173     hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
1174     hfetrack.SetRecTrack(track);
1175     if(!fPIDpreselect->IsSelected(&hfetrack)) {
1176       //printf("Did not pass AliHFEcuts::kPID\n");
1177       survived = kFALSE;
1178     }
1179     //else printf("Pass AliHFEcuts::kPID\n");
1180   }
1181
1182   return survived;
1183
1184 }
1185 //____________________________________________________________
1186 void AliAnalysisTaskHFE::MakeEventContainer(){
1187   //
1188   // Create the event container for the correction framework and link it
1189   // 1st bin: Vertex z-position
1190   // 2nd bin: V0AND decision (normalization to sigma_inel)
1191   // 3rd bin: Centrality class (for pp defined as 99.)
1192   //
1193   const Int_t kNvar = 3;  // number of variables on the grid:
1194   Int_t nBins[kNvar] = {120, 2, 11};
1195   Double_t binMin[kNvar] = {-30. , 0., 0.0};
1196   Double_t binMax[kNvar] = {30., 2., 11.0};
1197
1198   AliCFContainer *evCont = new AliCFContainer("eventContainer", "Container for events", AliHFEcuts::kNcutStepsEvent, kNvar, nBins);
1199
1200   Double_t *vertexBins = AliHFEtools::MakeLinearBinning(nBins[0], binMin[0], binMax[0]);
1201   Double_t *v0andBins = AliHFEtools::MakeLinearBinning(nBins[1], binMin[1], binMax[1]);
1202   Double_t *centralityBins = AliHFEtools::MakeLinearBinning(nBins[2], binMin[2], binMax[2]);
1203   evCont->SetBinLimits(0, vertexBins);
1204   evCont->SetBinLimits(1, v0andBins);
1205   evCont->SetBinLimits(2, centralityBins);
1206   delete[] vertexBins; delete[] v0andBins; delete[] centralityBins;
1207
1208   fCFM->SetEventContainer(evCont);
1209 }
1210
1211 //____________________________________________________________
1212 void AliAnalysisTaskHFE::MakeParticleContainer(){
1213   //
1214   // Create the particle container for the correction framework manager and
1215   // link it
1216   //
1217
1218   if(!fContainer) fContainer = new AliHFEcontainer("trackContainer");
1219   fVarManager->DefineVariables(fContainer);
1220
1221   // Create Correction Framework containers
1222   fContainer->CreateContainer("MCTrackCont", "Track Container filled with MC information", AliHFEcuts::kNcutStepsMCTrack);
1223   fContainer->CreateContainer("recTrackContReco", "Track Container filled with MC information", AliHFEcuts::kNcutStepsRecTrack + fPID->GetNumberOfPIDdetectors());
1224   fContainer->CreateContainer("recTrackContMC", "Track Container filled with MC information", AliHFEcuts::kNcutStepsRecTrack + fPID->GetNumberOfPIDdetectors());
1225
1226   fContainer->CreateContainer("hadronicBackground", "Container for Hadronic Background", 3);
1227   fContainer->CreateContainer("recTrackContDEReco", "Container for displaced electron analysis with Reco information", 1);
1228   fContainer->CreateContainer("recTrackContDEMC", "Container for displaced electron analysis with MC information", 1);
1229   fContainer->CreateCorrelationMatrix("correlationstepafterPID","THnSparse with correlations");
1230   fContainer->CreateCorrelationMatrix("correlationstepafterDE","THnSparse with correlations");
1231   if(!fVarManager->IsVariableDefined("centrality")) {
1232     //printf("Create the two other correlation maps\n");
1233     fContainer->CreateCorrelationMatrix("correlationstepbeforePID","THnSparse with correlations");
1234     fContainer->CreateCorrelationMatrix("correlationstepafterTOF","THnSparse with correlations");
1235   }
1236
1237   // Define the step names
1238   for(UInt_t istep = 0; istep < AliHFEcuts::kNcutStepsMCTrack; istep++){
1239     fContainer->SetStepTitle("MCTrackCont", AliHFEcuts::MCCutName(istep), istep);
1240   }
1241   for(UInt_t istep = 0; istep < AliHFEcuts::kNcutStepsRecTrack; istep++){
1242     fContainer->SetStepTitle("recTrackContReco", AliHFEcuts::RecoCutName(istep), istep);
1243     fContainer->SetStepTitle("recTrackContMC", AliHFEcuts::RecoCutName(istep), istep);
1244   }
1245   for(UInt_t ipid = 0; ipid < fPID->GetNumberOfPIDdetectors(); ipid++){
1246     fContainer->SetStepTitle("recTrackContReco", fPID->SortedDetectorName(ipid), AliHFEcuts::kNcutStepsRecTrack + ipid);
1247     fContainer->SetStepTitle("recTrackContMC", fPID->SortedDetectorName(ipid), AliHFEcuts::kNcutStepsRecTrack + ipid);
1248   }
1249 }
1250
1251 //____________________________________________________________
1252 void AliAnalysisTaskHFE::InitPIDperformanceQA(){
1253   // Add a histogram for Fake electrons
1254   const Int_t nDim=6;
1255   Int_t nBin[nDim] = {40, 8, 18, 2, 3, 3};
1256   //number of variables on the grid:pt,eta,phi,charge,
1257   const Double_t kPtbound[2] = {0.1, 20.};
1258   const Double_t kEtabound[2] = {-0.8, 0.8};
1259   const Double_t kPhibound[2] = {0., 2. * TMath::Pi()};
1260   const Double_t kChargebound[2] = {-1.1, 1.1};
1261   const Double_t kAddInf1bound[2] = {0., 3.};
1262   const Double_t kAddInf2bound[2] = {0., 3.};
1263   Double_t minima[nDim] = {kPtbound[0], kEtabound[0], kPhibound[0], kChargebound[0], kAddInf1bound[0], kAddInf2bound[0]};
1264   Double_t maxima[nDim] = {kPtbound[1], kEtabound[1], kPhibound[1], kChargebound[1], kAddInf1bound[1], kAddInf2bound[1]};
1265
1266   fQACollection->CreateTHnSparse("PIDperformance", "PID performance; pT [GeV/c]; theta [rad]; phi [rad]; charge; type (0 - not el, 1 - other el, 2 - HF el; flavor (0 - no, 1 - charm, 2 - bottom)", nDim, nBin, minima, maxima);
1267   fQACollection->CreateTHnSparse("SignalToBackgroundMC", "PID performance; pT [GeV/c]; theta [rad]; phi [rad]; charge; flavor (0 - no, 1 - charm, 2 - bottom); ITS Cluster (0 - no, 1 - first (and maybe second), 2 - second)", nDim, nBin, minima, maxima);
1268
1269   fQACollection->BinLogAxis("PIDperformance", 0);
1270   fQACollection->BinLogAxis("SignalToBackgroundMC", 0);
1271   fQACollection->Sumw2("PIDperformance");
1272   fQACollection->Sumw2("SignalToBackgroundMC");
1273 }
1274
1275 //____________________________________________________________
1276 void AliAnalysisTaskHFE::InitContaminationQA(){
1277   const Double_t kPtbound[2] = {0.1, 20.};
1278   Int_t iBin[1];
1279   iBin[0] = 44; // bins in pt
1280   fQACollection->CreateTH1F("hadronsBeforeIPcut", "Hadrons before IP cut", iBin[0], kPtbound[0], kPtbound[1], 1);
1281   fQACollection->CreateTH1F("hadronsAfterIPcut", "Hadrons after IP cut", iBin[0], kPtbound[0], kPtbound[1], 1);
1282   fQACollection->CreateTH1F("hadronsBeforeIPcutMC", "Hadrons before IP cut: MC p_{t}", iBin[0], kPtbound[0], kPtbound[1], 1);
1283   fQACollection->CreateTH1F("hadronsAfterIPcutMC", "Hadrons after IP cut: MC p_{t} ", iBin[0],kPtbound[0], kPtbound[1], 1);
1284 }
1285
1286 //____________________________________________________________
1287 void AliAnalysisTaskHFE::PrintStatus() const {
1288   //
1289   // Print Analysis status
1290   //
1291   printf("\n\tAnalysis Settings\n\t========================================\n\n");
1292   printf("\tSecondary Vertex finding: %s\n", GetPlugin(kSecVtx) ? "YES" : "NO");
1293   printf("\tPrimary Vertex resolution: %s\n", GetPlugin(kPriVtx) ? "YES" : "NO");
1294   printf("\tDisplaced electron analysis step: %s\n", GetPlugin(kDEstep) ? "YES" : "NO");
1295   printf("\tTagged Track Analysis: %s\n", GetPlugin(kTaggedTrackAnalysis) ? "YES" : "NO");
1296   printf("\n");
1297   printf("\tParticle Identification Detectors:\n");
1298   fPID->PrintStatus();
1299   printf("\n");
1300   printf("\tQA: \n");
1301   printf("\t\tPID: %s\n", IsQAOn(kPIDqa) ? "YES" :  "NO");
1302   printf("\t\tCUTS: %s\n", (fCuts != NULL && fCuts->IsQAOn()) ? "YES" : "NO");
1303   printf("\t\tMC: %s\n", IsQAOn(kMCqa) ? "YES" : "NO");
1304   printf("\n");
1305 }
1306
1307 //____________________________________________________________
1308 Bool_t AliAnalysisTaskHFE::FillProductionVertex(const AliVParticle * const track) const{
1309   //
1310   // Find the production vertex of the associated MC track
1311   //
1312   if(!fMCEvent) return kFALSE;
1313   const AliVParticle *mctrack = NULL;
1314   TString objectType = track->IsA()->GetName();
1315   if(objectType.CompareTo("AliESDtrack") == 0 || objectType.CompareTo("AliAODTrack") == 0){
1316     // Reconstructed track
1317     mctrack = fMCEvent->GetTrack(TMath::Abs(track->GetLabel()));
1318   } else {
1319     // MCParticle
1320     mctrack = track;
1321   }
1322
1323   if(!mctrack) return kFALSE;
1324
1325   Double_t xv = 0.0;
1326   Double_t yv = 0.0;
1327
1328   if(TString(mctrack->IsA()->GetName()).CompareTo("AliMCParticle") == 0){
1329     // case MCParticle
1330     const AliMCParticle *mcpart = dynamic_cast<const AliMCParticle *>(mctrack);
1331     if(mcpart){
1332       xv =  mcpart->Xv();
1333       yv =  mcpart->Yv();
1334     }
1335   } else {
1336     // case AODMCParticle
1337     const AliAODMCParticle *mcpart = dynamic_cast<const AliAODMCParticle *>(mctrack);
1338     if(mcpart){
1339       xv =  mcpart->Xv();
1340       yv =  mcpart->Yv();
1341     }
1342   }
1343
1344   //printf("xv %f, yv %f\n",xv,yv);
1345   fQACollection->Fill("radius", TMath::Abs(xv),TMath::Abs(yv));
1346
1347   return kTRUE;
1348
1349 }
1350 //__________________________________________
1351 void AliAnalysisTaskHFE::SwitchOnPlugin(Int_t plug){
1352   //
1353   // Switch on Plugin
1354   // Available:
1355   //  - Primary vertex studies
1356   //  - Secondary vertex Studies
1357   //  - Post Processing
1358   //
1359   switch(plug){
1360     case kPriVtx: SETBIT(fPlugins, plug); break;
1361     case kSecVtx: SETBIT(fPlugins, plug); break;
1362     case kIsElecBackGround: SETBIT(fPlugins, plug); break;
1363     case kPostProcess: SETBIT(fPlugins, plug); break;
1364     case kDEstep: SETBIT(fPlugins, plug); break;
1365     case kTaggedTrackAnalysis: SETBIT(fPlugins, plug); break;
1366     default: AliError("Unknown Plugin");
1367   };
1368 }
1369 //__________________________________________
1370 Bool_t AliAnalysisTaskHFE::ProcessCutStep(Int_t cutStep, AliVParticle *track){
1371   //
1372   // Check single track cuts for a given cut step
1373   // Fill the particle container
1374   //
1375   const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
1376   if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
1377   if(fVarManager->IsSignalTrack()) {
1378     fVarManager->FillContainer(fContainer, "recTrackContReco", cutStep, kFALSE);
1379     fVarManager->FillContainer(fContainer, "recTrackContMC", cutStep, kTRUE);
1380   }
1381   return kTRUE;
1382 }
1383 //___________________________________________________
1384 void AliAnalysisTaskHFE::ReadCentrality() {
1385   //
1386   // Recover the centrality of the event from ESD or AOD
1387   //
1388  if(IsAODanalysis()){
1389
1390    AliAODEvent *fAOD = dynamic_cast<AliAODEvent *>(fInputEvent);
1391    if(!fAOD){
1392      AliError("AOD Event required for AOD Analysis")
1393        return;
1394    }
1395    // Centrality
1396    //AliAODCentrality *aodCentrality = fAOD->GetCentrality();
1397    //Double_t fCentralityF = aodCentrality->GetCentralityPercentile("V0M");
1398    fCentralityF = 99.0; // Fake for the moment
1399
1400
1401  } else {
1402
1403    AliDebug(3, "Processing ESD Centrality");
1404    AliESDEvent *fESD = dynamic_cast<AliESDEvent *>(fInputEvent);
1405    if(!fESD){
1406      AliError("ESD Event required for ESD Analysis")
1407        return;
1408    }
1409    const char* type = fESD->GetBeamType();
1410
1411    if (strstr(type,"Pb-Pb")) {
1412
1413    // Centrality
1414    AliCentrality *esdCentrality = fESD->GetCentrality();
1415    Float_t fCentralityF_temp = esdCentrality->GetCentralityPercentile("V0M");
1416
1417    if( fCentralityF_temp >=  0. && fCentralityF_temp <   5.) fCentralityF =  0;
1418    else if ( fCentralityF_temp >=  5. && fCentralityF_temp <  10.) fCentralityF =  1;
1419    else if ( fCentralityF_temp >= 10. && fCentralityF_temp <  20.) fCentralityF = 2;
1420    else if ( fCentralityF_temp >= 20. && fCentralityF_temp <  30.) fCentralityF = 3;
1421    else if ( fCentralityF_temp >= 30. && fCentralityF_temp <  40.) fCentralityF = 4;
1422    else if ( fCentralityF_temp >= 40. && fCentralityF_temp <  50.) fCentralityF = 5;
1423    else if ( fCentralityF_temp >= 50. && fCentralityF_temp <  60.) fCentralityF = 6;
1424    else if ( fCentralityF_temp >= 60. && fCentralityF_temp <  70.) fCentralityF = 7;
1425    else if ( fCentralityF_temp >= 70. && fCentralityF_temp <  80.) fCentralityF = 8;
1426    else if ( fCentralityF_temp >= 80. && fCentralityF_temp <  90.) fCentralityF = 9;
1427    else if ( fCentralityF_temp >= 90. && fCentralityF_temp <=100.) fCentralityF = 10;
1428
1429    }
1430
1431
1432    if (strstr(type,"p-p")) {
1433      fCentralityF = 0;
1434      Int_t centralityF_temp = 0;
1435      const AliESDVertex *vtxESD = fESD->GetPrimaryVertexTracks();
1436      if(vtxESD && vtxESD->GetStatus()) centralityF_temp = vtxESD->GetNContributors();
1437
1438      //printf("centralityF_temp %d\n",centralityF_temp);
1439
1440      if( centralityF_temp <=  0) fCentralityF =  0;
1441      else if ( centralityF_temp >  0 && centralityF_temp <  2) fCentralityF = 1;
1442      else if ( centralityF_temp >=  2 && centralityF_temp <  3) fCentralityF = 2;
1443      else if ( centralityF_temp >= 3 && centralityF_temp <  4) fCentralityF = 3;
1444      else if ( centralityF_temp >= 4 && centralityF_temp <  5) fCentralityF = 4;
1445      else if ( centralityF_temp >= 5 && centralityF_temp <  10) fCentralityF = 5;
1446      else if ( centralityF_temp >= 10 && centralityF_temp <  20) fCentralityF = 6;
1447      else if ( centralityF_temp >= 20 && centralityF_temp <  30) fCentralityF = 7;
1448      else if ( centralityF_temp >= 30 && centralityF_temp <  40) fCentralityF = 8;
1449      else if ( centralityF_temp >= 40 && centralityF_temp <  50) fCentralityF = 9;
1450      else if ( centralityF_temp >= 50) fCentralityF = 10;
1451
1452
1453    }
1454
1455   //printf("centrality %f\n",fCentralityF);
1456
1457  }
1458
1459 }
1460 //___________________________________________________
1461 void AliAnalysisTaskHFE::RejectionPileUpVertexRangeEventCut() {
1462   //
1463   // Recover the centrality of the event from ESD or AOD
1464   //
1465  if(IsAODanalysis()){
1466
1467    AliAODEvent *fAOD = dynamic_cast<AliAODEvent *>(fInputEvent);
1468    if(!fAOD){
1469      AliError("AOD Event required for AOD Analysis")
1470        return;
1471    }
1472    // PileUp
1473    if(fRemovePileUp && fAOD->IsPileupFromSPD()) fIdentifiedAsPileUp = kTRUE;
1474    // Z vertex
1475    if(TMath::Abs(fAOD->GetPrimaryVertex()->GetZ()) > fCuts->GetVertexRange()) fIdentifiedAsOutInz = kTRUE;
1476    // Event Cut
1477    fPassTheEventCut = kTRUE;
1478    if(!fCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fAOD)) fPassTheEventCut = kFALSE;
1479
1480
1481  } else {
1482
1483    AliDebug(3, "Processing ESD Centrality");
1484    AliESDEvent *fESD = dynamic_cast<AliESDEvent *>(fInputEvent);
1485    if(!fESD){
1486      AliError("ESD Event required for ESD Analysis")
1487        return;
1488    }
1489    // PileUp
1490    fIdentifiedAsPileUp = kFALSE;
1491    if(fRemovePileUp && fESD->IsPileupFromSPD()) fIdentifiedAsPileUp = kTRUE;
1492    // Z vertex
1493    fIdentifiedAsOutInz = kFALSE;
1494    if(fESD->GetPrimaryVertexTracks()){
1495      if(TMath::Abs(fESD->GetPrimaryVertexTracks()->GetZ()) > fCuts->GetVertexRange()) fIdentifiedAsOutInz = kTRUE;
1496    }
1497    //Event Cut
1498    fPassTheEventCut = kTRUE;
1499    if(!fCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fESD)) fPassTheEventCut = kFALSE;
1500
1501
1502  }
1503
1504 }
1505