PWGCF/EBYE/Fluctuations/AliEbyEHigherMomentsTaskPID.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: Satyajit Jena.                                                 *
   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
  17 //=========================================================================//
  18 //                                                                         //
  19 //           Analysis Task for Net-Charge Higher Moment Analysis           //
  20 //              Author: Satyajit Jena || Nirbhay K. Behera                 //
  21 //                      sjena@cern.ch || nbehera@cern.ch                   //
  22 //                                                                         //
  23 //=========================================================================//
  24 #include "TChain.h"
  25 #include "TList.h"
  26 #include "TFile.h"
  27 #include "TTree.h"
  28 #include "TH1D.h"
  29 #include "TH2D.h"
  30 #include "TH3D.h"
  31 #include "THnSparse.h"
  32 #include "TCanvas.h"
  33
  34 #include "AliAnalysisTaskSE.h"
  35 #include "AliAnalysisManager.h"
  36 #include "AliAODEvent.h"
  37 #include "AliVTrack.h"
  38 #include "AliAODTrack.h"
  39 #include "AliAODInputHandler.h"
  40 #include "AliAODEvent.h"
  41 #include "AliStack.h"
  42 #include "AliAODMCHeader.h"
  43 #include "AliAODMCParticle.h"
  44 #include "AliMCEventHandler.h"
  45 #include "AliMCEvent.h"
  46 #include "AliStack.h"
  47 #include "AliGenHijingEventHeader.h"
  48 #include "AliGenEventHeader.h"
  49 #include "AliPID.h"
  50 #include "AliAODPid.h"
  51 #include "AliPIDResponse.h"
  52 #include "AliAODpidUtil.h"
  53 #include "AliPIDCombined.h"
  54 #include "AliHelperPID.h"
  55
  56 #include "AliEbyEHigherMomentsTaskPID.h"
  57
  58 using std::cout;
  59 using std::endl;
  60
  61 ClassImp(AliEbyEHigherMomentsTaskPID)
  62
  63
  64 //-----------------------------------------------------------------------
  65 AliEbyEHigherMomentsTaskPID::AliEbyEHigherMomentsTaskPID( const char *name )
  66 : AliAnalysisTaskSE( name ),
  67   fListOfHistos(0),
  68   fArrayMC(0),
  69   fPIDResponse(0),
  70   fParticleSpecies(AliPID::kProton),
  71   fAnalysisType(0),
  72   fCentralityEstimator("V0M"),
  73   fCentrality(0),
  74   fVxMax(3.),
  75   fVyMax(3.),
  76   fVzMax(10.),
  77   fPtLowerLimit(0.3),
  78   fPtHigherLimit(1.5),
  79   fNptBins(0),
  80   fBin(0),
  81   fEtaLowerLimit(-0.8),
  82   fEtaHigherLimit(0.8),
  83   fRapidityCut(0.5),
  84   fNSigmaCut(3.),
  85   fAODtrackCutBit(128),
  86   fHelperPID(0),
  87   fEventCounter(0),
  88   fHistVxVy(0),
  89   fTHnCentNplusNminusPid(0),
  90   fTHnCentNplusNminusPidTruth(0),
  91   fPtBinNplusNminusPid(0),
  92   fPtBinNplusNminusPidTruth(0)
  93 {
  94
  95   for ( Int_t i = 0; i < 4; i++) {
  96     fHistQA[i] = NULL;
  97   }
  98
  99   DefineOutput(1, TList::Class()); // Outpput...
 100 }
 101
 102 AliEbyEHigherMomentsTaskPID::~AliEbyEHigherMomentsTaskPID(){
 103
 104   if(fListOfHistos) delete fListOfHistos;
 105   if(fHelperPID) delete fHelperPID;
 106 }
 107
 108 //---------------------------------------------------------------------------------
 109 void AliEbyEHigherMomentsTaskPID::UserCreateOutputObjects() {
 110
 111   fListOfHistos = new TList();
 112   fListOfHistos->SetOwner(kTRUE);
 113
 114   fEventCounter = new TH1D("fEventCounter","EventCounter", 10, 0.5,10.5);
 115   fEventCounter->GetXaxis()->SetBinLabel(1,"Event Accesed");
 116   fEventCounter->GetXaxis()->SetBinLabel(2,"Within 0-80% centrality");
 117   fEventCounter->GetXaxis()->SetBinLabel(5,"Have a vertex");
 118   fEventCounter->GetXaxis()->SetBinLabel(6,"After vertex Cut");
 119   fEventCounter->GetXaxis()->SetBinLabel(7,"Event Analyzed");
 120   fListOfHistos->Add(fEventCounter);
 121
 122   //For QA-Histograms
 123   fHistQA[0] = new TH1D("fHistQAvz", "Histo Vz After Cut", 500, -25., 25.);
 124   fHistQA[1] = new TH1D("fHistQAPt","p_{T} distribution",600,0,6);
 125   fHistQA[2] = new TH1D("fHistQAEta","#eta distribution",240,-1.2,1.2);
 126   fHistQA[3] = new TH1D("fHistQAPhi","#phi distribution",340,0,6.8);
 127
 128   for(Int_t i = 0; i < 4; i++)
 129     {
 130       fListOfHistos->Add(fHistQA[i]);
 131     }
 132
 133   fHistVxVy = new TH2D("fHistVxVy","Vertex-x Vs Vertex-y", 200, -1., 1., 200, -1., 1.);
 134   fListOfHistos->Add(fHistVxVy);
 135
 136
 137   const Int_t nDim = 3;
 138   const Int_t nPid = 5;
 139   TString Species[nPid] = {"Electron","Muon","Pion","Kaon","Proton"};
 140
 141   Int_t fBins[nPid][nDim];
 142   Double_t fMin[nPid][nDim];
 143   Double_t fMax[nPid][nDim];
 144
 145   for( Int_t i = 0; i < nPid; i++ ){
 146     for( Int_t j = 0; j < nDim; j++ ){
 147       fBins[i][j] = 0;
 148       fMin[i][j] = 0.;
 149       fMax[i][j] = 0.;
 150     }
 151   }
 152
 153
 154   const Int_t dim = fNptBins*2 + 1;
 155   Int_t bin[nPid][dim];
 156   Double_t min[nPid][dim];
 157   Double_t max[nPid][dim];
 158
 159   bin[2][0] = 100; min[2][0] = -0.5; max[2][0] = 99.5;
 160   bin[3][0] = 100; min[3][0] = -0.5; max[3][0] = 99.5;
 161   bin[4][0] = 100; min[4][0] = -0.5; max[4][0] = 99.5;
 162
 163   for(Int_t i = 1; i < dim; i++){
 164
 165     bin[2][i] = 900;
 166     min[2][i] = -0.5;
 167     max[2][i] = 899.5;
 168
 169     bin[3][i] = 700;
 170     min[3][i] = -0.5;
 171     max[3][i] = 699.5;
 172
 173     bin[4][i] = 400;
 174     min[4][i] = -0.5;
 175     max[4][i] = 399.5;
 176
 177
 178   }
 179
 180
 181
 182
 183
 184   TString hname1, hname11;
 185   TString htitle1, axisTitle1, axisTitle2;
 186
 187   Int_t gPid = (Int_t)fParticleSpecies;
 188
 189   //Pion----
 190   fBins[2][0] = 100; fBins[2][1] = 1000; fBins[2][2] = 1000;
 191   fMin[2][0] = -0.5; fMin[2][1] = -0.5; fMin[2][2] = -0.5;
 192   fMax[2][0] = 99.5; fMax[2][1] = 999.5; fMax[2][2] = 999.5;
 193   //Kaon------
 194   fBins[3][0] = 100; fBins[3][1] = 700; fBins[3][2] = 700;
 195   fMin[3][0] = -0.5; fMin[3][1] = -0.5; fMin[3][2] = -0.5;
 196   fMax[3][0] = 99.5; fMax[3][1] = 699.5; fMax[3][2] = 699.5;
 197   //Proton-----
 198   fBins[4][0] = 100; fBins[4][1] = 400; fBins[4][2] = 400;
 199   fMin[4][0] = -0.5; fMin[4][1] = -0.5; fMin[4][2] = -0.5;
 200   fMax[4][0] = 99.5; fMax[4][1] = 399.5; fMax[4][2] = 399.5;
 201
 202   hname1 = "fCentNplusNminusPid"; hname1 +=gPid;
 203   htitle1 = Species[gPid]; htitle1 +=" And Neg-"; htitle1 +=Species[gPid];
 204   axisTitle1 ="Pos-"; axisTitle1 += Species[gPid];
 205   axisTitle2 ="Neg-"; axisTitle2 += Species[gPid];
 206
 207   fTHnCentNplusNminusPid = new THnSparseD(hname1.Data(),htitle1.Data(), nDim, fBins[gPid], fMin[gPid],fMax[gPid]);
 208
 209   fTHnCentNplusNminusPid->GetAxis(0)->SetTitle("Centrality");
 210   fTHnCentNplusNminusPid->GetAxis(1)->SetTitle(axisTitle1.Data());
 211   fTHnCentNplusNminusPid->GetAxis(2)->SetTitle(axisTitle2.Data());
 212   fListOfHistos->Add(fTHnCentNplusNminusPid);
 213
 214   TString hname2 = "fPtBinNplusNminusPid";
 215   TString htitle2 = "cent-nplus-nminus-ptbinwise";
 216   hname2 += gPid;
 217   htitle2 += gPid;
 218   fPtBinNplusNminusPid = new THnSparseI(hname2.Data(),htitle2.Data(), dim, bin[gPid], min[gPid], max[gPid]);
 219   fListOfHistos->Add(fPtBinNplusNminusPid);
 220
 221
 222
 223
 224   if( fAnalysisType == "MCAOD" ){
 225
 226     hname11 = "fCentNplusNminusPidTruth"; hname11 +=gPid;
 227     fTHnCentNplusNminusPidTruth = new THnSparseD(hname11.Data(),htitle1.Data(), nDim, fBins[gPid], fMin[gPid],fMax[gPid]);
 228
 229     fTHnCentNplusNminusPidTruth->GetAxis(0)->SetTitle("Centrality");
 230     fTHnCentNplusNminusPidTruth->GetAxis(1)->SetTitle(axisTitle1.Data());
 231     fTHnCentNplusNminusPidTruth->GetAxis(2)->SetTitle(axisTitle2.Data());
 232     fListOfHistos->Add(fTHnCentNplusNminusPidTruth);
 233
 234     TString hname22 = "fPtBinNplusNminusPidTruth";
 235     TString htitle22 = "cent-nplus-nminus-ptbinwise";
 236     hname22 += gPid;
 237     htitle22 += gPid;
 238
 239     fPtBinNplusNminusPidTruth = new THnSparseI(hname22.Data(),htitle22.Data(), dim, bin[gPid], min[gPid], max[gPid]);
 240     fListOfHistos->Add(fPtBinNplusNminusPidTruth);
 241
 242   }//fUsePid-------
 243
 244   PostData(1, fListOfHistos);
 245
 246
 247 }
 248
 249
 250 //----------------------------------------------------------------------------------
 251 void AliEbyEHigherMomentsTaskPID::UserExec( Option_t * ){
 252
 253   fEventCounter->Fill(1);
 254
 255   if(fAnalysisType == "AOD") {
 256
 257     doAODEvent();
 258
 259   }//AOD--analysis-----
 260
 261   else if(fAnalysisType == "MCAOD") {
 262
 263     doMCAODEvent();
 264
 265   }
 266
 267   else return;
 268
 269
 270
 271 }
 272
 273 //--------------------------------------------------------------------------------------
 274 void AliEbyEHigherMomentsTaskPID::doAODEvent(){
 275
 276   Double_t posPidSum = 0.;
 277   Double_t negPidSum = 0.;
 278
 279   const Int_t dim = fNptBins*2;
 280   Double_t PtCh[dim];
 281
 282   for(Int_t idx = 0; idx < dim; idx++){
 283     PtCh[idx] = 0.;
 284   }
 285
 286
 287   Int_t gPid = (Int_t)fParticleSpecies;
 288
 289   AliAnalysisManager* manager = AliAnalysisManager::GetAnalysisManager();
 290   if (!manager) {
 291     cout<<"ERROR: Analysis manager not found."<<endl;
 292     return;
 293   }
 294   //coneect to the inputHandler------------
 295   AliAODInputHandler* inputHandler = dynamic_cast<AliAODInputHandler*> (manager->GetInputEventHandler());
 296   if (!inputHandler) {
 297     cout<<"ERROR: Input handler not found."<<endl;
 298     return;
 299   }
 300
 301   AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
 302
 303   if (!fAOD)
 304     {
 305       cout<< "ERROR: AOD not found " <<endl;
 306       return;
 307     }
 308
 309   fPIDResponse = inputHandler->GetPIDResponse();
 310
 311
 312   if (!fPIDResponse){
 313     AliFatal("This Task needs the PID response attached to the inputHandler");
 314     return;
 315   }
 316
 317   if(!ProperVertex(fAOD)) return;
 318
 319   AliAODHeader *aodHeader = dynamic_cast<AliAODHeader*>(fAOD->GetHeader());
 320   if(!aodHeader) AliFatal("Not a standard AOD");
 321
 322   fCentrality = (Int_t)aodHeader->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data());
 323
 324   if(fCentrality < 0 || fCentrality >= 81) return;
 325
 326   fEventCounter->Fill(2);
 327
 328   Int_t nTracks = fAOD->GetNumberOfTracks();
 329
 330   for(Int_t i = 0; i < nTracks; i++) {
 331
 332     AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(fAOD->GetTrack(i));
 333
 334     if(!aodTrack) {
 335       AliError(Form("ERROR: Could not retrieve AODtrack %d",i));
 336       continue;
 337     }
 338
 339     if(!AcceptTrack(aodTrack) ) continue;
 340
 341
 342     fBin = GetPtBin(aodTrack->Pt());
 343
 344     Short_t gCharge = aodTrack->Charge();
 345
 346
 347     Double_t rap =  aodTrack->Y(AliPID::ParticleMass(fParticleSpecies));
 348     if( fabs(rap ) > fRapidityCut ) continue;//Rapidity cut
 349
 350     Int_t PID = fHelperPID->GetParticleSpecies((AliVTrack*)aodTrack, kFALSE);
 351
 352     if( (PID+2) == gPid ){
 353
 354       if(gCharge > 0){
 355         posPidSum += 1.;
 356         PtCh[fBin] += 1;
 357       }
 358       if(gCharge < 0){
 359         negPidSum += 1.;
 360         PtCh[fNptBins+fBin] += 1;
 361       }
 362
 363     }//chek the PID
 364
 365
 366   }//--------- Track Loop to select with filterbit
 367
 368
 369   Double_t fContainerPid[3] = { static_cast<Double_t>(fCentrality), posPidSum, negPidSum};
 370   fTHnCentNplusNminusPid->Fill(fContainerPid);
 371
 372   //cout << fCentrality <<" "<< posPidSum <<" " << negPidSum << endl;
 373
 374   Double_t ptContainer[dim+1];
 375
 376   ptContainer[0] = fCentrality;
 377
 378   for(Int_t i = 1; i <= dim; i++){
 379     ptContainer[i] = PtCh[i-1];
 380     //cout << PtCh[i-1] <<" ,";
 381   }
 382   //cout << endl;
 383
 384   fPtBinNplusNminusPid->Fill(ptContainer);
 385
 386
 387   fEventCounter->Fill(7);
 388   PostData(1, fListOfHistos);
 389   return;
 390
 391 }
 392 //--------------------------------------------------------------------------------------
 393 //--------------------------------------------------------------------------------------
 394 void AliEbyEHigherMomentsTaskPID::doMCAODEvent(){
 395
 396
 397   //---------
 398   Double_t posPidSumMCRec = 0.;
 399   Double_t negPidSumMCRec = 0.;
 400
 401   Double_t posPidSumMCTruth = 0.;
 402   Double_t negPidSumMCTruth = 0.;
 403
 404   const Int_t dim = fNptBins*2;
 405   Double_t PtCh[dim];
 406   Double_t ptChMC[dim];
 407
 408   for(Int_t idx = 0; idx < dim; idx++){
 409     PtCh[idx] = 0;
 410     ptChMC[idx] = 0;
 411   }
 412
 413
 414   Int_t gPid = (Int_t)fParticleSpecies;
 415
 416   Int_t gPdgCode = AliPID::ParticleCode(fParticleSpecies);
 417
 418   //Connect to Anlaysis manager------
 419   AliAnalysisManager* manager = AliAnalysisManager::GetAnalysisManager();
 420   if (!manager) {
 421     cout<<"ERROR: Analysis manager not found."<<endl;
 422     return;
 423   }
 424
 425   AliAODInputHandler* inputHandler = dynamic_cast<AliAODInputHandler*> (manager->GetInputEventHandler());
 426   if (!inputHandler) {
 427     cout<<"ERROR: Input handler not found."<<endl;
 428     return;
 429   }
 430
 431   //AOD event------
 432   AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
 433
 434   if (!fAOD)
 435     {
 436       cout<< "ERROR: AOD not found " <<endl;
 437       return;
 438     }
 439
 440   fPIDResponse =inputHandler->GetPIDResponse();
 441
 442
 443   if(!fPIDResponse){
 444     AliFatal("This Task needs the PID response attached to the inputHandler");
 445     return;
 446   }
 447
 448   // -- Get MC header
 449   // ------------------------------------------------------------------
 450
 451   fArrayMC = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
 452   if (!fArrayMC) {
 453     AliFatal("Error: MC particles branch not found!\n");
 454     return;
 455   }
 456
 457   AliAODMCHeader *mcHdr=NULL;
 458   mcHdr=(AliAODMCHeader*)fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName());
 459   if(!mcHdr) {
 460     printf("MC header branch not found!\n");
 461     return;
 462   }
 463
 464   if(!ProperVertex(fAOD)) return;
 465
 466   AliAODHeader *aodHeader = dynamic_cast<AliAODHeader*>(fAOD->GetHeader());
 467   if(!aodHeader) AliFatal("Not a standard AOD");
 468
 469   fCentrality = (Int_t)aodHeader->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data());
 470
 471   if( fCentrality < 0 || fCentrality >= 81) return;
 472
 473   fEventCounter->Fill(2);
 474
 475   Int_t nTracks = fAOD->GetNumberOfTracks();
 476
 477   for(Int_t i = 0; i < nTracks; i++) {
 478
 479     AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(fAOD->GetTrack(i));
 480
 481     if(!aodTrack) {
 482       AliError(Form("ERROR: Could not retrieve AODtrack %d",i));
 483       continue;
 484     }
 485
 486     if(!AcceptTrack(aodTrack) ) continue;
 487
 488     fBin = GetPtBin(aodTrack->Pt());
 489
 490     //cout << "Pt Bin " << fBin << endl;
 491
 492     Short_t gCharge = aodTrack->Charge();
 493
 494     Double_t rap =  aodTrack->Y(AliPID::ParticleMass(fParticleSpecies));
 495     if( fabs(rap ) > fRapidityCut ) continue;//Rapidity cut
 496
 497     Int_t PID = fHelperPID->GetParticleSpecies((AliVTrack*)aodTrack, kFALSE);
 498
 499     if( (PID+2) == gPid  ){
 500       if(gCharge > 0){
 501         posPidSumMCRec += 1;
 502         PtCh[fBin] += 1;
 503       }
 504       if( gCharge < 0 ){
 505         negPidSumMCRec += 1.;
 506         PtCh[fNptBins+fBin] += 1;
 507       }
 508     }//nSigmaCut-----
 509
 510   }//--------- Track Loop to select with filterbit
 511
 512   //===============================================================================
 513   //---------------------MC Truth--------------------------------------------------
 514   //===============================================================================
 515
 516   Int_t nMCTrack = fArrayMC->GetEntriesFast();
 517
 518   for (Int_t iMC = 0; iMC < nMCTrack; iMC++) {
 519
 520     AliAODMCParticle *partMC = (AliAODMCParticle*) fArrayMC->At(iMC);
 521
 522     if(!partMC){
 523       AliError(Form("ERROR: Could not retrieve AODMCtrack %d",iMC));
 524       continue;
 525     }
 526
 527     if(partMC->Charge() == 0) continue;
 528     if(!partMC->IsPhysicalPrimary()) continue;
 529
 530     if (partMC->Eta() < fEtaLowerLimit || partMC->Eta() > fEtaHigherLimit) continue;
 531     if (partMC->Pt() < fPtLowerLimit ||  partMC->Pt() > fPtHigherLimit) continue;
 532
 533     Short_t gCharge = partMC->Charge();
 534     Int_t mcbin = GetPtBin(partMC->Pt());
 535
 536     Double_t rap =  partMC->Y();
 537     if( fabs(rap ) > fRapidityCut ) continue;//Rapidity cut
 538     if(TMath::Abs(partMC->GetPdgCode()) != gPdgCode) continue;
 539
 540     if(gCharge > 0){
 541       posPidSumMCTruth += 1.;
 542       ptChMC[mcbin] += 1;
 543     }
 544     if(gCharge < 0){
 545       negPidSumMCTruth += 1.;
 546       ptChMC[fNptBins+mcbin] += 1;
 547     }
 548
 549
 550   }//MC-Truth Track loop--
 551
 552
 553   //cout << fCentrality << " " << positiveSumMCRec << " " << negativeSumMCRec <<endl;
 554   //cout <<"   " << positiveSumMCTruth << " " << negativeSumMCTruth << endl;
 555   //cout <<fCentrality<<"   " << posPidSumMCRec << " " << negPidSumMCRec << endl;
 556   //cout <<fCentrality<<"   " << posPidSumMCTruth << " " << negPidSumMCTruth << endl;
 557   //cout <<"---------------------------------" << endl;
 558
 559   Double_t fContainerPid[3] = { static_cast<Double_t>(fCentrality), posPidSumMCRec, negPidSumMCRec};//Reco.
 560   Double_t fContainerPidTruth[3] = { static_cast<Double_t>(fCentrality), posPidSumMCTruth, negPidSumMCTruth};
 561
 562   fTHnCentNplusNminusPid->Fill(fContainerPid);//Fill the rec. pid tracks
 563   fTHnCentNplusNminusPidTruth->Fill(fContainerPidTruth);//MC-Truth pid
 564
 565   Double_t ptContainer[dim+1];
 566   Double_t ptContainerMC[dim+1];
 567   ptContainer[0] = fCentrality;
 568   ptContainerMC[0] = fCentrality;
 569
 570   for(Int_t i = 1; i <= dim; i++){
 571     ptContainer[i] = PtCh[i-1];
 572     ptContainerMC[i] = ptChMC[i-1];
 573     //cout <<" "<< PtCh[i-1]<<endl;
 574     //cout<< " Rec=" << PtCh[i-1];
 575   }
 576
 577   //cout << endl;
 578
 579   fPtBinNplusNminusPid->Fill(ptContainer);
 580   fPtBinNplusNminusPidTruth->Fill(ptContainerMC);
 581
 582   fEventCounter->Fill(7);
 583   PostData(1, fListOfHistos);
 584   return;
 585
 586 }
 587
 588 //---------------------------------------------------------------------------------------
 589 Bool_t AliEbyEHigherMomentsTaskPID::ProperVertex(AliAODEvent *fAOD) const{
 590
 591   Bool_t IsVtx = kFALSE;
 592
 593   const AliAODVertex *vertex = fAOD->GetPrimaryVertex();
 594
 595   if(vertex) {
 596     Double32_t fCov[6];
 597     vertex->GetCovarianceMatrix(fCov);
 598     if(vertex->GetNContributors() > 0) {
 599       if(fCov[5] != 0) {
 600
 601         Double_t lvx = vertex->GetX();
 602         Double_t lvy = vertex->GetY();
 603         Double_t lvz = vertex->GetZ();
 604
 605         fEventCounter->Fill(5);
 606
 607         if(TMath::Abs(lvx) < fVxMax) {
 608           if(TMath::Abs(lvy) < fVyMax) {
 609             if(TMath::Abs(lvz) < fVzMax) {
 610
 611               fEventCounter->Fill(6);
 612               fHistQA[0]->Fill(lvz);
 613               fHistVxVy->Fill(lvx,lvy);
 614               IsVtx = kTRUE;
 615
 616             }//Z-Vertex cut---
 617           }//Y-vertex cut--
 618         }//X-vertex cut---
 619       }//Covariance------
 620     }//Contributors check---
 621   }//If vertex-----------
 622
 623   AliCentrality *centrality = fAOD->GetCentrality();
 624   if (centrality->GetQuality() != 0) IsVtx = kFALSE;
 625
 626   return IsVtx;
 627 }
 628 //---------------------------------------------------------------------------------------
 629
 630 //---------------------------------------------------------------------------------------
 631 Bool_t AliEbyEHigherMomentsTaskPID::AcceptTrack(AliAODTrack* track) const{
 632
 633   if(!track) return kFALSE;
 634   if( track->Charge() == 0 ) return kFALSE;
 635
 636   Double_t pt = track->Pt();
 637   Double_t eta = track->Eta();
 638   if(!track->TestFilterBit(fAODtrackCutBit)) return kFALSE;
 639   if( pt < fPtLowerLimit || pt > fPtHigherLimit ) return kFALSE;
 640   if( eta < fEtaLowerLimit || eta > fEtaHigherLimit ) return kFALSE;
 641
 642   fHistQA[1]->Fill(pt);
 643   fHistQA[2]->Fill(eta);
 644   fHistQA[3]->Fill(track->Phi());
 645
 646   return kTRUE;
 647 }
 648 //------------------------------------------------------------------------
 649
 650 //------------------------------------------------------------------------
 651 Int_t AliEbyEHigherMomentsTaskPID::GetPtBin(Double_t pt){
 652
 653   Int_t bin = 0;
 654
 655   Double_t BinSize = (fPtHigherLimit - fPtLowerLimit)/fNptBins;
 656
 657   for(Int_t iBin = 0; iBin < fNptBins; iBin++){
 658
 659     Double_t xLow = fPtLowerLimit + iBin*BinSize;
 660     Double_t xHigh = fPtLowerLimit + (iBin+1)*BinSize;
 661
 662     if( pt >= xLow && pt < xHigh){
 663       bin = iBin;
 664       //cout << "Pt Bin #" << bin <<"pt=" << pt << endl;
 665       break;
 666     }
 667
 668   }//for
 669
 670
 671   return bin;
 672 }
 673 //------------------------------------------------------------------------
 674 //------------------------------------------------------------------------
 675 void AliEbyEHigherMomentsTaskPID::Terminate( Option_t * ){
 676
 677   Info("AliEbyEHigherMomentTask"," Task Successfully finished");
 678
 679 }
 680
 681 //------------------------------------------------------------------------