PWGCF/EBYE/NetParticle/AliAnalysisNetParticleHelper.cxx

   1 //-*- Mode: C++ -*-
   2
   3 #include "TMath.h"
   4 #include "TAxis.h"
   5 #include "TSystem.h"
   6 #include "TFile.h"
   7 #include "TPRegexp.h"
   8
   9 #include "AliStack.h"
  10 #include "AliMCEvent.h"
  11 #include "AliMCParticle.h"
  12 #include "AliESDtrackCuts.h"
  13 #include "AliESDInputHandler.h"
  14 #include "AliESDpid.h"
  15 #include "AliCentrality.h"
  16 #include "AliTracker.h"
  17
  18 #include "AliAnalysisNetParticleHelper.h"
  19
  20 using namespace std;
  21
  22 // Task for NetParticle checks
  23 // Author: Jochen Thaeder <jochen@thaeder.de>
  24
  25 ClassImp(AliAnalysisNetParticleHelper)
  26
  27 /*
  28  * ---------------------------------------------------------------------------------
  29  *                            Constructor / Destructor
  30  * ---------------------------------------------------------------------------------
  31  */
  32
  33 //________________________________________________________________________
  34 AliAnalysisNetParticleHelper::AliAnalysisNetParticleHelper() :
  35   fESDHandler(NULL),
  36   fPIDResponse(NULL),
  37   fESD(NULL),
  38   fMCEvent(NULL),
  39   fStack(NULL),
  40
  41   fCentralityBin(-1),
  42   fCentralityPercentile(-1.),
  43
  44   fCentralityBinMax(7),
  45   fVertexZMax(10.),
  46   fRapidityMax(0.5),
  47   fMinTrackLengthMC(70.),
  48   fNSigmaMaxCdd(3.),
  49   fNSigmaMaxCzz(3.),
  50
  51   fParticleSpecies(AliPID::kProton),
  52   fControlParticleCode(3122),
  53   fControlParticleIsNeutral(kTRUE),
  54   fControlParticleName("Lambda"),
  55
  56   fNSigmaMaxTPC(2.5),
  57   fNSigmaMaxTOF(2.5),
  58   fMinPtForTOFRequired(0.8),
  59
  60   fHEventStat0(NULL),
  61   fHEventStat1(NULL),
  62   fHEventStatMax(6),
  63
  64   fHTriggerStat(NULL),
  65   fNTriggers(5),
  66
  67   fHCentralityStat(NULL),
  68   fNCentralityBins(10),
  69
  70   fEtaCorrFunc(NULL),
  71   fCorr0(NULL),
  72   fCorr1(NULL) {
  73   // Constructor
  74
  75   AliLog::SetClassDebugLevel("AliAnalysisNetParticleHelper",10);
  76 }
  77
  78 //________________________________________________________________________
  79 AliAnalysisNetParticleHelper::~AliAnalysisNetParticleHelper() {
  80   // Destructor
  81
  82   for (Int_t jj = 0; jj < 2; ++jj) {
  83     if (fCorr0[jj]) delete[] fCorr0[jj];
  84     if (fCorr1[jj]) delete[] fCorr1[jj];
  85   }
  86   if (fCorr0) delete[] fCorr0;
  87   if (fCorr1) delete[] fCorr1;
  88
  89   return;
  90 }
  91
  92 /*
  93  * ---------------------------------------------------------------------------------
  94  *                                 Public Methods
  95  * ---------------------------------------------------------------------------------
  96  */
  97
  98 //________________________________________________________________________
  99 Int_t AliAnalysisNetParticleHelper::Initialize(Bool_t isMC) {
 100   // -- Initialize helper
 101
 102   Int_t iResult = 0;
 103
 104   // -- Setup event cut statistics
 105   InitializeEventStats();
 106
 107   // -- Setup trigger statistics
 108   InitializeTriggerStats();
 109
 110   // -- Setup centrality statistics
 111   InitializeCentralityStats();
 112
 113   // -- Load Eta correction function
 114   iResult = InitializeEtaCorrection(isMC);
 115
 116   // -- Load Eta correction function
 117   iResult = InitializeTrackbyTrackCorrection();
 118
 119   return iResult;
 120 }
 121
 122 //________________________________________________________________________
 123 Int_t AliAnalysisNetParticleHelper::SetupEvent(AliESDInputHandler *esdHandler, AliMCEvent *mcEvent) {
 124   // -- Setup Event
 125
 126   // -- Get ESD objects
 127   fESDHandler  = esdHandler;
 128   fPIDResponse = esdHandler->GetPIDResponse();
 129   fESD         = fESDHandler->GetEvent();
 130
 131   // -- Get MC objects
 132   fMCEvent     = mcEvent;
 133   if (fMCEvent)
 134     fStack     = fMCEvent->Stack();
 135
 136   // -- Get event centrality
 137   // >  0-5|5-10|10-20|20-30|30-40|40-50|50-60|60-70|70-80|80-90 --> 10 bins
 138   // >   0   1     2     3     4     5     6     7     8     9
 139
 140   AliCentrality *centrality = fESD->GetCentrality();
 141
 142   Int_t centBin = centrality->GetCentralityClass10("V0M");
 143   if (centBin == 0)
 144     fCentralityBin = centrality->GetCentralityClass5("V0M");
 145   else if (centBin == 10 || centBin == -1.)
 146     fCentralityBin = -1;
 147   else if (centBin > 0 && centBin < fNCentralityBins)
 148     fCentralityBin = centBin + 1;
 149   else
 150     fCentralityBin = -2;
 151
 152   // -- Stay within the max centrality bin
 153   if (fCentralityBin >= fCentralityBinMax)
 154     fCentralityBin = -2;
 155
 156   fCentralityPercentile = centrality->GetCentralityPercentile("V0M");
 157
 158   // -- Update TPC pid with eta correction
 159   UpdateEtaCorrectedTPCPid();
 160
 161   return 0;
 162 }
 163
 164 /*
 165  * ---------------------------------------------------------------------------------
 166  *                         Event / Trigger Statistics
 167  * ---------------------------------------------------------------------------------
 168  */
 169
 170 //________________________________________________________________________
 171 Bool_t AliAnalysisNetParticleHelper::IsEventTriggered() {
 172   // -- Check if Event is triggered and fill Trigger Histogram
 173
 174   Bool_t *aTriggerFired = new Bool_t[fNTriggers];
 175   for (Int_t ii = 0; ii < fNTriggers; ++ii)
 176     aTriggerFired[ii] = kFALSE;
 177
 178   if ((fESDHandler->IsEventSelected() & AliVEvent::kMB))          aTriggerFired[0] = kTRUE;
 179   if ((fESDHandler->IsEventSelected() & AliVEvent::kCentral))     aTriggerFired[1] = kTRUE;
 180   if ((fESDHandler->IsEventSelected() & AliVEvent::kSemiCentral)) aTriggerFired[2] = kTRUE;
 181   if ((fESDHandler->IsEventSelected() & AliVEvent::kEMCEJE))      aTriggerFired[3] = kTRUE;
 182   if ((fESDHandler->IsEventSelected() & AliVEvent::kEMCEGA))      aTriggerFired[4] = kTRUE;
 183
 184   Bool_t isTriggered = kFALSE;
 185
 186   for (Int_t ii=0; ii<fNTriggers; ++ii) {
 187     if(aTriggerFired[ii]) {
 188       isTriggered = kTRUE;
 189       fHTriggerStat->Fill(ii);
 190     }
 191   }
 192
 193   delete[] aTriggerFired;
 194
 195   return isTriggered;
 196 }
 197
 198 //________________________________________________________________________
 199 Bool_t AliAnalysisNetParticleHelper::IsEventRejected() {
 200   // -- Evaluate event statistics histograms
 201
 202   Int_t *aEventCuts = new Int_t[fHEventStatMax];
 203   // set aEventCuts[ii] to 1 in case of reject
 204
 205   for (Int_t ii=0;ii<fHEventStatMax; ++ii)
 206     aEventCuts[ii] = 0;
 207
 208   Int_t iCut = 0;
 209
 210   // -- 0 - Before Physics Selection
 211   aEventCuts[iCut] = 0;
 212
 213   // -- 1 - No Trigger fired
 214   ++iCut;
 215   if (!IsEventTriggered())
 216     aEventCuts[iCut] = 1;
 217
 218   // -- 2 - No Vertex
 219   ++iCut;
 220   const AliESDVertex* vtxESD = fESD->GetPrimaryVertexTracks();
 221   if (!vtxESD)
 222     aEventCuts[iCut] = 1;
 223
 224   // -- 3 - Vertex z outside cut window
 225   ++iCut;
 226   if (vtxESD){
 227     if(TMath::Abs(vtxESD->GetZv()) > fVertexZMax)
 228       aEventCuts[iCut] = 1;
 229   }
 230   else
 231     aEventCuts[iCut] = 1;
 232
 233   // -- 4 - Centrality = -1  (no centrality or not hadronic)
 234   ++iCut;
 235   if(fCentralityBin == -1.)
 236     aEventCuts[iCut] = 1;
 237
 238   // -- 5 - Centrality < fCentralityMax
 239   ++iCut;
 240   if(fCentralityBin == -2.)
 241     aEventCuts[iCut] = 1;
 242
 243   // -- Fill statistics / reject event
 244   Bool_t isRejected = FillEventStats(aEventCuts);
 245
 246   // -- Cleanup
 247   delete[] aEventCuts;
 248
 249   return isRejected;
 250 }
 251
 252 /*
 253  * ---------------------------------------------------------------------------------
 254  *                         Accept Particle Methods - private
 255  * ---------------------------------------------------------------------------------
 256  */
 257
 258 //________________________________________________________________________
 259 Bool_t AliAnalysisNetParticleHelper::IsParticleAcceptedBasicCharged(TParticle *particle, Int_t idxMC) {
 260   // -- Check if MC particle is accepted for basic parameters
 261
 262   if (!particle)
 263     return kFALSE;
 264
 265   // -- check if PDF code exists
 266   if (!particle->GetPDG())
 267     return kFALSE;
 268
 269   // -- check if charged
 270   if (particle->GetPDG()->Charge() == 0.0)
 271     return kFALSE;
 272
 273   // -- check if physical primary
 274   if(!fStack->IsPhysicalPrimary(idxMC))
 275     return kFALSE;
 276
 277   return kTRUE;
 278 }
 279 //________________________________________________________________________
 280 Bool_t AliAnalysisNetParticleHelper::IsParticleAcceptedBasicNeutral(TParticle *particle, Int_t idxMC) {
 281   // -- Check if MC particle is accepted for basic parameters
 282
 283   if (!particle)
 284     return kFALSE;
 285
 286   // -- check if PDF code exists
 287   if (!particle->GetPDG())
 288     return kFALSE;
 289
 290   // -- check if neutral
 291   if (particle->GetPDG()->Charge() != 0.0)
 292     return kFALSE;
 293
 294   // -- check if physical primary
 295   if(!fStack->IsPhysicalPrimary(idxMC))
 296     return kFALSE;
 297
 298   return kTRUE;
 299 }
 300
 301 //________________________________________________________________________
 302 Bool_t AliAnalysisNetParticleHelper::IsParticleAcceptedRapidity(TParticle *particle, Double_t &yP) {
 303   // -- Check if particle is accepted
 304   // > in rapidity
 305   // > return 0 if not accepted
 306
 307   Double_t mP = AliPID::ParticleMass(fParticleSpecies);
 308
 309   // -- Calculate rapidities and kinematics
 310   Double_t p  = particle->P();
 311   Double_t pz = particle->Pz();
 312
 313   Double_t eP = TMath::Sqrt(p*p + mP*mP);
 314   yP          = 0.5 * TMath::Log((eP + pz) / (eP - pz));
 315
 316   // -- Check Rapidity window
 317   if (TMath::Abs(yP) > fRapidityMax)
 318     return kFALSE;
 319
 320   return kTRUE;
 321 }
 322
 323 //_____________________________________________________________________________
 324 Bool_t AliAnalysisNetParticleHelper::IsParticleFindable(Int_t label) {
 325   // -- Check if MC particle is findable tracks
 326
 327   AliMCParticle *mcParticle = static_cast<AliMCParticle*>(fMCEvent->GetTrack(label));
 328   if(!mcParticle)
 329     return kFALSE;
 330
 331   Int_t counter;
 332   Float_t tpcTrackLength = mcParticle->GetTPCTrackLength(AliTracker::GetBz(), 0.05, counter, 3.0);
 333
 334   return (tpcTrackLength > fMinTrackLengthMC);
 335 }
 336
 337 /*
 338  * ---------------------------------------------------------------------------------
 339  *                            Accept Track Methods - public
 340  * ---------------------------------------------------------------------------------
 341  */
 342
 343 //________________________________________________________________________
 344 Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedBasicCharged(AliESDtrack *track) {
 345   // -- Check if track is accepted
 346   // > for basic parameters
 347
 348   if (!track)
 349     return kFALSE;
 350
 351   if (track->Charge() == 0)
 352     return kFALSE;
 353
 354   // -- Get momentum for dEdx
 355   if (!track->GetInnerParam())
 356     return kFALSE;
 357
 358   return kTRUE;
 359 }
 360
 361 //________________________________________________________________________
 362 Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedRapidity(AliESDtrack *track, Double_t &yP) {
 363   // -- Check if track is accepted
 364   // > in rapidity
 365   // > return 0 if not accepted
 366
 367   Double_t mP = AliPID::ParticleMass(fParticleSpecies);
 368
 369   // -- Calculate rapidities and kinematics
 370   Double_t pvec[3];
 371   track->GetPxPyPz(pvec);
 372
 373   Double_t p  = track->GetP();
 374   Double_t eP = TMath::Sqrt(p*p + mP*mP);
 375            yP = 0.5 * TMath::Log((eP + pvec[2]) / (eP - pvec[2]));
 376
 377   // -- Check Rapidity window
 378   if (TMath::Abs(yP) > fRapidityMax)
 379     return kFALSE;
 380
 381   return kTRUE;
 382 }
 383
 384 //________________________________________________________________________
 385 Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedDCA(AliESDtrack *track) {
 386   // -- Check if track is accepted
 387   // > for DCA, if both SPD layers have hits
 388
 389   Bool_t isAccepted = kTRUE;
 390
 391   // -- Get nHits SPD
 392   if (track->HasPointOnITSLayer(0) && track->HasPointOnITSLayer(1)) {
 393
 394     // -- Get DCA nSigmas
 395     Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
 396     track->GetImpactParameters(dca,cov);
 397
 398     Float_t nSigmaCdd = (cov[0] != 0.) ? dca[0]/TMath::Sqrt(cov[0]) : -9.99;
 399     Float_t nSigmaCzz = (cov[2] != 0.) ? dca[1]/TMath::Sqrt(cov[2]) : -9.99;
 400
 401     if (fNSigmaMaxCdd != 0.) {
 402       if (TMath::Abs(nSigmaCdd) > fNSigmaMaxCdd)
 403         isAccepted = kFALSE;
 404     }
 405
 406     if (fNSigmaMaxCzz != 0.) {
 407       if (TMath::Abs(nSigmaCzz) > fNSigmaMaxCzz)
 408         isAccepted = kFALSE;
 409     }
 410   }
 411
 412   return isAccepted;
 413 }
 414
 415 //________________________________________________________________________
 416 Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedPID(AliESDtrack *track, Double_t* pid) {
 417   // -- Check if track is accepted
 418   // > provides TPC and TOF nSigmas to the argument
 419
 420   Bool_t isAcceptedTPC = kFALSE;
 421   Bool_t isAcceptedTOF = kTRUE;
 422
 423   // -- Get PID
 424   pid[0] = fPIDResponse->NumberOfSigmasTPC(track, fParticleSpecies);
 425   pid[1] = fPIDResponse->NumberOfSigmasTOF(track, fParticleSpecies);
 426
 427   // -- Check PID with TPC
 428   if (TMath::Abs(pid[0]) < fNSigmaMaxTPC)
 429     isAcceptedTPC = kTRUE;
 430
 431   // -- Check PID with TOF
 432   if (isAcceptedTPC) {
 433
 434     // Has TOF PID availible
 435     if (track->GetStatus() & AliVTrack::kTOFpid) {
 436       if (TMath::Abs(pid[1]) < fNSigmaMaxTOF)
 437         isAcceptedTOF = kTRUE;
 438       else
 439         isAcceptedTOF = kFALSE;
 440     }
 441
 442     // Has no TOF PID availible
 443     else {
 444       if (track->Pt() > fMinPtForTOFRequired)
 445         isAcceptedTOF = kFALSE;
 446       else
 447         isAcceptedTOF = kTRUE;
 448     }
 449   } // if (isAcceptedTPC) {
 450
 451   return (isAcceptedTPC && isAcceptedTOF);
 452 }
 453
 454 /*
 455  * ---------------------------------------------------------------------------------
 456  *                         Helper Methods
 457  * ---------------------------------------------------------------------------------
 458  */
 459
 460 //________________________________________________________________________
 461 void AliAnalysisNetParticleHelper::UpdateEtaCorrectedTPCPid() {
 462   // -- Update eta corrected TPC pid
 463
 464   if (!fEtaCorrFunc)
 465     return;
 466
 467   for (Int_t idxTrack = 0; idxTrack < fESD->GetNumberOfTracks(); ++idxTrack) {
 468     AliESDtrack *track = fESD->GetTrack(idxTrack);
 469
 470     // -- Check if charged track is accepted for basic parameters
 471     if (!IsTrackAcceptedBasicCharged(track))
 472       continue;
 473
 474     Double_t newTPCSignal = track->GetTPCsignal() / fEtaCorrFunc->Eval(track->Eta());
 475     track->SetTPCsignal(newTPCSignal, track->GetTPCsignalSigma(), track->GetTPCsignalN());
 476   }
 477 }
 478
 479 //________________________________________________________________________
 480 Double_t AliAnalysisNetParticleHelper::GetTrackbyTrackCorrectionFactor(Double_t *aTrack,  Int_t flag) {
 481   // -- Get efficiency correctionf of particle dependent on (eta, phi, pt, centrality)
 482
 483   Int_t idxPart = (aTrack[2] < 0) ? 0 : 1;
 484   THnSparseF* corrMatrix = (flag == 0) ? fCorr0[idxPart][fCentralityBin] : fCorr1[idxPart][fCentralityBin];
 485
 486   Double_t dimBin[3] = {aTrack[3], aTrack[4], aTrack[1]}; // eta, phi, pt
 487
 488   Double_t corr = corrMatrix->GetBinContent(corrMatrix->GetBin(dimBin));
 489   if (corr == 0.) {
 490     AliError(Form("Should not happen : bin content = 0. >> eta: %.2f | phi : %.2f | pt : %.2f | cent %d",
 491                   aTrack[3], aTrack[4], aTrack[1], fCentralityBin));
 492     corr = 1.;
 493   }
 494
 495   return corr;
 496 }
 497
 498 //________________________________________________________________________
 499 void AliAnalysisNetParticleHelper::BinLogAxis(const THnSparseF *h, Int_t axisNumber) {
 500   // -- Method for the correct logarithmic binning of histograms
 501
 502   TAxis *axis = h->GetAxis(axisNumber);
 503   Int_t  bins = axis->GetNbins();
 504
 505   Double_t from  = axis->GetXmin();
 506   Double_t to    = axis->GetXmax();
 507   Double_t *newBins = new Double_t[bins + 1];
 508
 509   newBins[0] = from;
 510   Double_t factor = pow(to/from, 1./bins);
 511
 512   for (int i = 1; i <= bins; i++) {
 513    newBins[i] = factor * newBins[i-1];
 514   }
 515   axis->Set(bins, newBins);
 516   delete [] newBins;
 517 }
 518
 519 ///////////////////////////////////////////////////////////////////////////////////
 520
 521 /*
 522  * ---------------------------------------------------------------------------------
 523  *                           Initialize - Private
 524  * ---------------------------------------------------------------------------------
 525  */
 526
 527 //________________________________________________________________________
 528 void AliAnalysisNetParticleHelper::InitializeEventStats() {
 529   // -- Initialize event statistics histograms
 530
 531   const Char_t* aEventNames[]      = {"All", "IsTriggered", "HasVertex", "Vz<Vz_{Max}", "Centrality [0,90]%"};
 532   const Char_t* aCentralityMaxNames[] = {"5", "10", "20", "30", "40", "50", "60", "70", "80", "90", "100"};
 533
 534   fHEventStat0 = new TH1F("fHEventStat0","Event cut statistics 0;Event Cuts;Events", fHEventStatMax,-0.5,fHEventStatMax-0.5);
 535   fHEventStat1 = new TH1F("fHEventStat1","Event cut statistics 1;Event Cuts;Events", fHEventStatMax,-0.5,fHEventStatMax-0.5);
 536
 537   for ( Int_t ii=0; ii < fHEventStatMax-1; ii++ ) {
 538     fHEventStat0->GetXaxis()->SetBinLabel(ii+1, aEventNames[ii]);
 539     fHEventStat1->GetXaxis()->SetBinLabel(ii+1, aEventNames[ii]);
 540   }
 541   fHEventStat0->GetXaxis()->SetBinLabel(fHEventStatMax, Form("Centrality [0-%s]%%", aCentralityMaxNames[fCentralityBinMax-1]));
 542   fHEventStat1->GetXaxis()->SetBinLabel(fHEventStatMax, Form("Centrality [0-%s]%%", aCentralityMaxNames[fCentralityBinMax-1]));
 543 }
 544
 545 //________________________________________________________________________
 546 void AliAnalysisNetParticleHelper::InitializeTriggerStats() {
 547   // -- Initialize trigger statistics histograms
 548
 549   const Char_t* aTriggerNames[] = { "kMB", "kCentral", "kSemiCentral", "kEMCEJE", "kEMCEGA" };
 550
 551   fHTriggerStat = new TH1F("fHTriggerStat","Trigger statistics;Trigger;Events", fNTriggers,-0.5,fNTriggers-0.5);
 552
 553   for ( Int_t ii=0; ii < fNTriggers; ii++ )
 554     fHTriggerStat->GetXaxis()->SetBinLabel(ii+1, aTriggerNames[ii]);
 555 }
 556
 557 //________________________________________________________________________
 558 void AliAnalysisNetParticleHelper::InitializeCentralityStats() {
 559   // -- Initialize trigger statistics histograms
 560
 561   const Char_t* aCentralityNames[] = {"0-5%", "5-10%", "10-20%", "20-30%", "30-40%", "40-50%",
 562                                       "50-60%", "60-70%", "70-80%", "80-90%", "90-100%"};
 563
 564   fHCentralityStat = new TH1F("fHCentralityStat","Centrality statistics;Centrality Bins;Events",
 565                               fNCentralityBins,-0.5,fNCentralityBins-0.5);
 566
 567   for ( Int_t ii=0; ii < fNCentralityBins; ii++ )
 568     fHCentralityStat->GetXaxis()->SetBinLabel(ii+1, aCentralityNames[ii]);
 569 }
 570
 571 //________________________________________________________________________
 572 Int_t AliAnalysisNetParticleHelper::InitializeEtaCorrection(Bool_t isMC) {
 573   // -- Initialize eta correction maps for TPC pid
 574
 575   TFile fileEtaCorrMaps("${ALICE_ROOT}/PWGDQ/dielectron/files/EtaCorrMaps.root");
 576   if (!fileEtaCorrMaps.IsOpen())
 577     return -1;
 578
 579   TList *keyList = fileEtaCorrMaps.GetListOfKeys();
 580
 581   TString sList("");
 582   sList = (isMC) ? "LHC11a10" :  "LHC10h.pass2";
 583
 584   for (Int_t idx = 0; idx < keyList->GetEntries(); ++idx) {
 585     TString keyName = keyList->At(idx)->GetName();
 586     TPRegexp reg(keyName);
 587     if (reg.MatchB(sList)){
 588       AliInfo(Form("Using Eta Correction Function: %s",keyName.Data()));
 589       fEtaCorrFunc = static_cast<TF1*>(fileEtaCorrMaps.Get(keyName.Data()));
 590       return 0;
 591     }
 592   }
 593
 594   return -2;
 595 }
 596
 597 //________________________________________________________________________
 598 Int_t AliAnalysisNetParticleHelper::InitializeTrackbyTrackCorrection() {
 599   // -- Initialize track by track correction matrices
 600
 601   AliInfo("TODO ... get the correct name from particle"); // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 602
 603   TFile* corrFile  = TFile::Open("${ALICE_ROOT}/PWGCF/EBYE/NetParticle/eff/effectiveCorrectionProton.root");
 604   if (!corrFile) {
 605     AliError("Track-by-track correction file can not be opened!");
 606     return -1;
 607   }
 608
 609   // -- correction - not cross section corrected
 610   fCorr0    = new THnSparseF**[2];
 611   fCorr0[0] = new THnSparseF*[fCentralityBinMax];
 612   fCorr0[1] = new THnSparseF*[fCentralityBinMax];
 613
 614   for (Int_t idxCent = 0; idxCent < fCentralityBinMax; ++idxCent) {
 615     THnSparseF *sp0 = static_cast<THnSparseF*>(corrFile->Get(Form("pbar_Corr0_Cent_%d", idxCent)));
 616     THnSparseF *sp1 = static_cast<THnSparseF*>(corrFile->Get(Form("p_Corr0_Cent_%d", idxCent)));
 617
 618     if (!sp0 || !sp1) {
 619       AliError(Form("Effective correction objects 0 - idxCent %d can not be retrieved!",idxCent));
 620       return -1;
 621     }
 622
 623     fCorr0[0][idxCent] = static_cast<THnSparseF*>(sp0->Clone());
 624     fCorr0[1][idxCent] = static_cast<THnSparseF*>(sp1->Clone());
 625   }
 626
 627   // -- correction - ross section corrected
 628   fCorr1    = new THnSparseF**[2];
 629   fCorr1[0] = new THnSparseF*[fCentralityBinMax];
 630   fCorr1[1] = new THnSparseF*[fCentralityBinMax];
 631
 632   for (Int_t idxCent = 0; idxCent < fCentralityBinMax; ++idxCent) {
 633     THnSparseF *sp0 = static_cast<THnSparseF*>(corrFile->Get(Form("pbar_Corr1_Cent_%d", idxCent)));
 634     THnSparseF *sp1 = static_cast<THnSparseF*>(corrFile->Get(Form("p_Corr1_Cent_%d", idxCent)));
 635
 636     if (!sp0 || !sp1) {
 637       AliError(Form("Effective correction objects 1 - idxCent %d can not be retrieved!",idxCent));
 638       return -1;
 639     }
 640
 641     fCorr1[0][idxCent] = static_cast<THnSparseF*>(sp0->Clone());
 642     fCorr1[1][idxCent] = static_cast<THnSparseF*>(sp1->Clone());
 643   }
 644
 645   return 0;
 646 }
 647
 648 /*
 649  * ---------------------------------------------------------------------------------
 650  *                         Event / Trigger Statistics - private
 651  * ---------------------------------------------------------------------------------
 652  */
 653
 654 //________________________________________________________________________
 655 Bool_t AliAnalysisNetParticleHelper::FillEventStats(Int_t *aEventCuts) {
 656   // -- Fill event / centrality statistics
 657
 658   Bool_t isRejected = kFALSE;
 659
 660   // -- Fill event statistics
 661   for (Int_t idx = 0; idx < fHEventStatMax ; ++idx) {
 662     if (aEventCuts[idx])
 663       isRejected = kTRUE;
 664     else
 665       fHEventStat0->Fill(idx);
 666   }
 667
 668   for (Int_t idx = 0; idx < fHEventStatMax; ++idx) {
 669     if (aEventCuts[idx])
 670       break;
 671     fHEventStat1->Fill(idx);
 672   }
 673
 674   // -- Fill centrality statistics of accepted events
 675   if (!isRejected)
 676     fHCentralityStat->Fill(fCentralityBin);
 677
 678   return isRejected;
 679 }