PWGPP/TPC/AliTPCPIDEtaQA.cxx

   1 #include "TChain.h"
   2 #include "TTree.h"
   3 #include "TF1.h"
   4 #include "TAxis.h"
   5 #include "TH1F.h"
   6 #include "THnSparse.h"
   7
   8 #include "AliMCParticle.h"
   9 //#include "AliStack.h"
  10
  11 #include "AliAnalysisTask.h"
  12 #include "AliAnalysisManager.h"
  13
  14 #include "AliESDEvent.h"
  15 #include "AliMCEvent.h"
  16 #include "AliESDInputHandler.h"
  17 #include "AliInputEventHandler.h"
  18
  19 #include "AliVVertex.h"
  20 #include "AliAnalysisFilter.h"
  21 #include "AliPID.h"
  22 #include "AliPIDResponse.h"
  23 #include "AliTPCPIDResponse.h"
  24
  25 #include "AliTPCPIDEtaQA.h"
  26
  27 /*
  28 This task determines the eta dependence of the TPC signal.
  29 For this purpose, only tracks fulfilling some standard quality cuts are taken into account.
  30 The obtained data can be used to derive the functional behaviour of the eta dependence.
  31 Such a function can be plugged into this task to correct for the eta dependence and to see
  32 if there is then really no eta dependence left.
  33
  34 Class written by Benjamin Hess.
  35 Contact: bhess@cern.ch
  36 */
  37
  38 ClassImp(AliTPCPIDEtaQA)
  39
  40 //________________________________________________________________________
  41 AliTPCPIDEtaQA::AliTPCPIDEtaQA()
  42   : AliTPCPIDBase()
  43   , fPtThresholdForPhiCut(2.0)
  44   , fPhiCutSecondBranchLow(0x0)
  45   , fPhiCutSecondBranchHigh(0x0)
  46   , fhPIDdataAll(0x0)
  47   //OLD clusterQA, fhNumClustersPhiPrimePtBeforeCut(0x0)
  48   //OLD clusterQA, fhNumClustersPhiPrimePtAfterCut(0x0)
  49   , fhPhiPrimeCutEfficiency(0x0)
  50   , fOutputContainer(0x0)
  51 {
  52   // default Constructor
  53
  54   // Question: Is this the right place to initialize these functions?
  55   // Will it work on proof? i.e. will they be streamed to the workers?
  56   // Also one should add getters and setters
  57   //TODO fPhiCutLow  = new TF1("StandardPhiCutLow",  "0.1/x/x+pi/18.0-0.025", 0, 50);
  58   //TODO fPhiCutHigh = new TF1("StandardPhiCutHigh", "0.12/x+pi/18.0+0.035", 0, 50);
  59
  60   //TODO NEW fPhiCutLow  = new TF1("StandardPhiCutLow",  "0.072/x+pi/18.0-0.035", 0, 50);
  61   //TODO NEW fPhiCutHigh = new TF1("StandardPhiCutHigh", "0.07/x/x+0.1/x+pi/18.0+0.035", 0, 50);
  62
  63   fPhiCutSecondBranchLow  = new TF1("StandardPhiCutSecondBranchLow",  "NewStandardPhiCutLow - 2.*pi/18.", 0, 50);
  64   fPhiCutSecondBranchHigh = new TF1("StandardPhiCutSecondBranchHigh", "0.07/x+pi/18.0+0.125 - 2.*pi/18.", 0, 50);
  65 }
  66
  67 //________________________________________________________________________
  68 AliTPCPIDEtaQA::AliTPCPIDEtaQA(const char *name)
  69   : AliTPCPIDBase(name)
  70   , fPtThresholdForPhiCut(2.0)
  71   , fPhiCutSecondBranchLow(0x0)
  72   , fPhiCutSecondBranchHigh(0x0)
  73   , fhPIDdataAll(0x0)
  74   //OLD clusterQA, fhNumClustersPhiPrimePtBeforeCut(0x0)
  75   //OLD clusterQA, fhNumClustersPhiPrimePtAfterCut(0x0)
  76   , fhPhiPrimeCutEfficiency(0x0)
  77   , fOutputContainer(0x0)
  78 {
  79   // Constructor
  80
  81   // Question: Is this the right place to initialize these functions?
  82   // Will it work on proof? i.e. will they be streamed to the workers?
  83   // Also one should add getters and setters
  84   //TODO fPhiCutLow  = new TF1("StandardPhiCutLow",  "0.1/x/x+pi/18.0-0.025", 0, 50);
  85   //TODO fPhiCutHigh = new TF1("StandardPhiCutHigh", "0.12/x+pi/18.0+0.035", 0, 50);
  86
  87   //TODO NEW fPhiCutLow  = new TF1("StandardPhiCutLow",  "0.072/x+pi/18.0-0.035", 0, 50);
  88   //TODO NEW fPhiCutHigh = new TF1("StandardPhiCutHigh", "0.07/x/x+0.1/x+pi/18.0+0.035", 0, 50);
  89
  90   fPhiCutSecondBranchLow  = new TF1("StandardPhiCutSecondBranchLow",  "StandardPhiCutLow - 2.*pi/18.", 0, 50);
  91   fPhiCutSecondBranchHigh = new TF1("StandardPhiCutSecondBranchHigh", "0.07/x+pi/18.0+0.125 - 2.*pi/18.", 0, 50);
  92
  93   // Define input and output slots here
  94   // Input slot #0 works with a TChain
  95   DefineInput(0, TChain::Class());
  96   // Output slot #1 writes into a TObjArray container
  97   DefineOutput(1, TObjArray::Class());
  98 }
  99
 100
 101 //________________________________________________________________________
 102 AliTPCPIDEtaQA::~AliTPCPIDEtaQA()
 103 {
 104   // dtor
 105
 106   delete fOutputContainer;
 107   fOutputContainer = 0;
 108
 109   delete fPhiCutSecondBranchLow;
 110   fPhiCutSecondBranchLow = 0;
 111
 112   delete fPhiCutSecondBranchHigh;
 113   fPhiCutSecondBranchHigh = 0;
 114 }
 115
 116
 117 //________________________________________________________________________
 118 void AliTPCPIDEtaQA::UserCreateOutputObjects()
 119 {
 120   // Create histograms
 121   // Called once
 122
 123   AliTPCPIDBase::UserCreateOutputObjects();
 124
 125   OpenFile(1);
 126
 127   fOutputContainer = new TObjArray(1);
 128   fOutputContainer->SetName(GetName()) ;
 129   fOutputContainer->SetOwner(kTRUE);
 130
 131   const Int_t nEta = TMath::Nint(40 * fEtaCut);
 132
 133   const Int_t nPtBins = 68;
 134   Double_t binsPt[nPtBins+1] = {0. ,  0.05, 0.1,  0.15, 0.2,  0.25, 0.3,  0.35, 0.4,  0.45,
 135            0.5,  0.55, 0.6,  0.65, 0.7,  0.75, 0.8,  0.85, 0.9,  0.95,
 136            1.0,  1.1 , 1.2,  1.3 , 1.4,  1.5 , 1.6,  1.7 , 1.8,  1.9 ,
 137            2.0,  2.2 , 2.4,  2.6 , 2.8,  3.0 , 3.2,  3.4 , 3.6,  3.8 ,
 138            4.0,  4.5 , 5.0,  5.5 , 6.0,  6.5 , 7.0,  8.0 , 9.0,  10.0,
 139            11.0, 12.0, 13.0, 14.0, 15.0, 16.0, 18.0, 20.0, 22.0, 24.0,
 140            26.0, 28.0, 30.0, 32.0, 34.0, 36.0, 40.0, 45.0, 50.0 };
 141
 142
 143   const Int_t nBins = 6;
 144   // MC PID, SelectSpecies, P(TPC_inner), multiplicity, deltaPrimeSpecies, eta
 145
 146   Int_t bins[nBins] =
 147     {  9,   4, nPtBins,    40, 201, nEta };
 148   Double_t xmin[nBins] =
 149     {  0., 0.,      0.,    0., 0.5, -fEtaCut};
 150   Double_t xmax[nBins] =
 151     {  9., 4.,    50.0, 20000, 2.0, fEtaCut };
 152
 153
 154   fhPIDdataAll = new THnSparseI("hPIDdataAll","", nBins, bins, xmin, xmax);
 155   SetUpHist(fhPIDdataAll, binsPt);
 156   fOutputContainer->Add(fhPIDdataAll);
 157
 158   /*OLD clusterQA
 159   const Int_t nBinsQA = 3;
 160   // pT, phiPrime, Ncl
 161
 162   Int_t binsQA[nBinsQA] =
 163     {  nPtBins, 50, 90 };
 164   Double_t xminQA[nBinsQA] =
 165     {  0., 0., 70.};
 166   Double_t xmaxQA[nBinsQA] =
 167     {  50.0, TMath::Pi() / 9., 160 };
 168
 169   fhNumClustersPhiPrimePtBeforeCut = new THnSparseI("hNumClustersPhiPrimePtBeforeCut", "", nBinsQA, binsQA, xminQA, xmaxQA);
 170   fhNumClustersPhiPrimePtBeforeCut->SetBinEdges(0, binsPt);
 171   fhNumClustersPhiPrimePtBeforeCut->GetAxis(0)->SetTitle("p_{T} (GeV/c)");
 172   fhNumClustersPhiPrimePtBeforeCut->GetAxis(1)->SetTitle("#phi'");
 173   fhNumClustersPhiPrimePtBeforeCut->GetAxis(2)->SetTitle("Ncl");
 174   fOutputContainer->Add(fhNumClustersPhiPrimePtBeforeCut);
 175
 176   fhNumClustersPhiPrimePtAfterCut = new THnSparseI("hNumClustersPhiPrimePtAfterCut", "", nBinsQA, binsQA, xminQA, xmaxQA);
 177   fhNumClustersPhiPrimePtAfterCut->SetBinEdges(0, binsPt);
 178   fhNumClustersPhiPrimePtAfterCut->GetAxis(0)->SetTitle("p_{T} (GeV/c)");
 179   fhNumClustersPhiPrimePtAfterCut->GetAxis(1)->SetTitle("#phi'");
 180   fhNumClustersPhiPrimePtAfterCut->GetAxis(2)->SetTitle("Ncl");
 181   fOutputContainer->Add(fhNumClustersPhiPrimePtAfterCut);
 182
 183   fhPhiPrimeCutEfficiency = new TH1F("hPhiPrimeCutEfficiency", "Efficiency of #phi' cut;p_{T} (GeV/c);Cut efficiency", nPtBins, binsPt);
 184   fOutputContainer->Add(fhPhiPrimeCutEfficiency);
 185   */
 186
 187   PostData(1, fOutputContainer);
 188 }
 189
 190
 191 //________________________________________________________________________
 192 void AliTPCPIDEtaQA::UserExec(Option_t *)
 193 {
 194   // Main loop
 195   // Called for each event
 196
 197   fESD = dynamic_cast<AliESDEvent*>(InputEvent());
 198   if (!fESD) {
 199     Printf("ERROR: fESD not available");
 200     return;
 201   }
 202
 203   fMC = dynamic_cast<AliMCEvent*>(MCEvent());
 204
 205   if (!fPIDResponse || !fV0KineCuts)
 206     return;
 207
 208   if (!GetVertexIsOk(fESD))
 209     return;
 210
 211   const AliVVertex* primaryVertex = fESD->GetPrimaryVertexTracks();
 212   if (!primaryVertex)
 213     return;
 214
 215   if(primaryVertex->GetNContributors() <= 0)
 216     return;
 217
 218   Double_t magField = fESD->GetMagneticField();
 219
 220   // Fill V0 arrays with V0s
 221   FillV0PIDlist();
 222
 223   Int_t multiplicity = fESD->GetNumberOfESDTracks();
 224
 225   // Track loop to fill a Train spectrum
 226   for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
 227     AliESDtrack* track = fESD->GetTrack(iTracks);
 228     if (!track) {
 229       Printf("ERROR: Could not receive track %d", iTracks);
 230       continue;
 231     }
 232
 233     if (TMath::Abs(track->Eta()) >= fEtaCut)  continue;
 234
 235     //if (TMath::Abs(track->Eta()) < 0.6 || TMath::Abs(track->Eta()) > 0.8)  continue;
 236
 237     // Do not distinguish positively and negatively charged V0's
 238     Char_t v0tagAllCharges = TMath::Abs(GetV0tag(iTracks));
 239     if (v0tagAllCharges == -99) {
 240       AliError(Form("Problem with V0 tag list (requested V0 track for track %d from %d (list status %d))!", iTracks, fESD->GetNumberOfTracks(),
 241                     fV0tags != 0x0));
 242       continue;
 243     }
 244
 245     Bool_t isV0el = v0tagAllCharges == 14;
 246     Bool_t isV0pi = v0tagAllCharges == 15;
 247     Bool_t isV0pr = v0tagAllCharges == 16;
 248     Bool_t isV0 = isV0el || isV0pi || isV0pr;
 249
 250     // Apply track cut. Accept track, if from V0 or if accepted by cut.
 251     // Do not apply the track cuts to V0s, since the track cuts are usually
 252     // cuts on primaries, what will throw away almost all V0s.
 253     if(!isV0 && (fTrackFilter && !fTrackFilter->IsSelected(track)))
 254       continue;
 255
 256     if (GetUseTPCCutMIGeo()) {
 257       // If cut on geometry is active, don't cut on number of clusters, since such a cut is already included.
 258       if (!isV0) {
 259         if (!TPCCutMIGeo(track, InputEvent()))
 260           continue;
 261       }
 262       else {
 263         // One should not cut on geometry for V0's since they have different topology. (Loosely) cut on num of clusters instead.
 264         if (track->GetTPCsignalN() < 60)
 265           continue;
 266       }
 267     }
 268     else {
 269       // If cut on geometry is not active, always cut on num clusters
 270       if (track->GetTPCsignalN() < 60)
 271         continue;
 272     }
 273
 274     Double_t pT = track->Pt();
 275     //Double_t phiPrime = GetPhiPrime(track->Phi(), magField, track->Charge());
 276
 277     //OLD clusterQA Double_t entryQA[3] = { pT, phiPrime, track->GetTPCsignalN() };
 278     //OLD clusterQA fhNumClustersPhiPrimePtBeforeCut->Fill(entryQA);
 279
 280     // Apply PhiPrimeCut only on high-pt tracks, in this case: Tracks with pt >= fPtThresholdForPhiCut
 281     if(fUsePhiCut && pT >= fPtThresholdForPhiCut) {
 282
 283       if(fUsePhiCut) {
 284       if (!PhiPrimeCut(track, magField))
 285         continue; // reject track
 286     }
 287       //Use cut for second branch?
 288       if(!PhiPrimeCut(track, magField))
 289       //if(!PhiPrimeCut(track, magField)
 290       //   || (phiPrime < fPhiCutSecondBranchHigh->Eval(pT)  && phiPrime > fPhiCutSecondBranchLow->Eval(pT)))
 291               continue; // reject track
 292     }
 293
 294     //OLD clusterQA fhNumClustersPhiPrimePtAfterCut->Fill(entryQA);
 295
 296     Int_t label = track->GetLabel();
 297
 298     AliMCParticle* mcTrack = 0x0;
 299
 300     if (fMC) {
 301       if (label < 0)
 302         continue; // If MC is available, reject tracks with negative ESD label
 303       mcTrack = dynamic_cast<AliMCParticle*>(fMC->GetTrack(TMath::Abs(label)));
 304       if (!mcTrack) {
 305         Printf("ERROR: Could not receive mcTrack with label %d for track %d", label, iTracks);
 306         continue;
 307       }
 308
 309       /*// Only accept MC primaries
 310       if (!fMC->Stack()->IsPhysicalPrimary(TMath::Abs(label))) {
 311         continue;
 312       }*/
 313     }
 314
 315     // Momentum
 316     Double_t pTPC = track->GetTPCmomentum();
 317
 318     // Eta
 319     Float_t eta = track->Eta();
 320
 321     // TPC signal
 322     Double_t dEdxTPC = track->GetTPCsignal();
 323
 324     Double_t dEdxExpectedEl = fPIDResponse->GetTPCResponse().GetExpectedSignal(track, AliPID::kElectron, AliTPCPIDResponse::kdEdxDefault,
 325                                                                               fPIDResponse->UseTPCEtaCorrection(),
 326                                                                               fPIDResponse->UseTPCMultiplicityCorrection());
 327     Double_t dEdxExpectedKa = fPIDResponse->GetTPCResponse().GetExpectedSignal(track, AliPID::kKaon, AliTPCPIDResponse::kdEdxDefault,
 328                                                                               fPIDResponse->UseTPCEtaCorrection(),
 329                                                                               fPIDResponse->UseTPCMultiplicityCorrection());
 330     Double_t dEdxExpectedPi = fPIDResponse->GetTPCResponse().GetExpectedSignal(track, AliPID::kPion, AliTPCPIDResponse::kdEdxDefault,
 331                                                                               fPIDResponse->UseTPCEtaCorrection(),
 332                                                                               fPIDResponse->UseTPCMultiplicityCorrection());
 333     Double_t dEdxExpectedPr = fPIDResponse->GetTPCResponse().GetExpectedSignal(track, AliPID::kProton, AliTPCPIDResponse::kdEdxDefault,
 334                                                                               fPIDResponse->UseTPCEtaCorrection(),
 335                                                                               fPIDResponse->UseTPCMultiplicityCorrection());
 336
 337     Double_t deltaPrimeEl = (dEdxExpectedEl > 0) ? (dEdxTPC / dEdxExpectedEl) : (-1);
 338     Double_t deltaPrimeKa = (dEdxExpectedKa > 0) ? (dEdxTPC / dEdxExpectedKa) : (-1);
 339     Double_t deltaPrimePi = (dEdxExpectedPi > 0) ? (dEdxTPC / dEdxExpectedPi) : (-1);
 340     Double_t deltaPrimePr = (dEdxExpectedPr > 0) ? (dEdxTPC / dEdxExpectedPr) : (-1);
 341
 342     Int_t binMC = -1;
 343     if (fMC) {
 344       Int_t pdg = mcTrack->PdgCode();
 345
 346       if (TMath::Abs(pdg) == 211) {//Pion
 347         // If V0, only accept if correctly identified
 348         if (isV0)   {
 349           if (isV0pi)
 350             binMC = 7;
 351           else
 352             continue;
 353         }
 354         else
 355           binMC = 3;
 356       }
 357       else if (TMath::Abs(pdg) == 321) {//Kaon
 358         // No kaons from V0 => Reject
 359         if (isV0)
 360           continue;
 361         else
 362           binMC = 1;
 363       }
 364       else if (TMath::Abs(pdg) == 2212) {//Proton
 365         // If V0, only accept if correctly identified
 366         if (isV0)   {
 367           if (isV0pr)
 368             binMC = 8;
 369           else
 370             continue;
 371         }
 372         else
 373           binMC = 4;
 374       }
 375       else if (TMath::Abs(pdg) == 11) {//Electron
 376         // If V0, only accept if correctly identified
 377         if (isV0)   {
 378           if (isV0el)
 379             binMC = 6;
 380           else
 381             continue;
 382         }
 383         else
 384           binMC = 0;
 385       }
 386       else if (TMath::Abs(pdg) == 13) {//Muon
 387         // No muons from V0 => Reject
 388         if (isV0)
 389           continue;
 390         else
 391           binMC = 2;
 392       }
 393       else
 394         continue; // Reject all other particles
 395     }
 396
 397     Bool_t isKaon = kFALSE;
 398     Bool_t isPion = kFALSE;
 399     Bool_t isProton = kFALSE;
 400     Bool_t isElectron = kFALSE;
 401     Bool_t isV0plusTOFel = kFALSE;
 402
 403
 404     if (!fMC && !isV0) {
 405       UInt_t status = track->GetStatus();
 406       Bool_t hasTOFout  = status&AliESDtrack::kTOFout;
 407       Bool_t hasTOFtime = status&AliESDtrack::kTIME;
 408       Bool_t hasTOF     = kFALSE;
 409       if (hasTOFout && hasTOFtime)
 410         hasTOF = kTRUE;
 411       Float_t length = track->GetIntegratedLength();
 412       // Length check only for primaries, not for V0's!
 413       if (length < 350 && !isV0)
 414         hasTOF = kFALSE;
 415
 416
 417       // Select Kaons, Pions and Protons in 3 sigma band (TOF and TPC). Below some momentum threshold, only use TPC cut,
 418       // since TOF efficiency is really bad.
 419       // In case of ambiguity, add entries for corresponding species
 420       if ((pTPC <= 0.25 && TMath::Abs(fPIDResponse->NumberOfSigmasTPC(track, AliPID::kKaon)) < 6.0) ||
 421           (pTPC > 0.25 && pTPC <= 0.3 && TMath::Abs(fPIDResponse->NumberOfSigmasTPC(track, AliPID::kKaon)) < 4.0) ||
 422           (pTPC > 0.3 && pTPC <= 0.35 && TMath::Abs(fPIDResponse->NumberOfSigmasTPC(track, AliPID::kKaon)) < 3.0) ||
 423           (pTPC > 0.35 && TMath::Abs(fPIDResponse->NumberOfSigmasTPC(track, AliPID::kKaon)) < 3.0 && hasTOF &&
 424            TMath::Abs(fPIDResponse->NumberOfSigmasTOF(track, AliPID::kKaon)) < 3.0)) {
 425           isKaon = kTRUE;
 426       }
 427       if ((dEdxTPC >= 50. / TMath::Power(pTPC, 1.3)) && // Pattern recognition instead of TPC cut to be ~independent of old TPC expected dEdx
 428          ((pTPC < 0.6) ||
 429          (pTPC >= 0.6 && hasTOF && TMath::Abs(fPIDResponse->NumberOfSigmasTOF(track, AliPID::kProton)) < 3.0))) {
 430         isProton = kTRUE;
 431       }
 432       /*
 433       if ((TMath::Abs(fPIDResponse->NumberOfSigmasTPC(track, AliPID::kProton)) < (pTPC < 0.3 ? 8.0 : 5.0)) &&
 434           ((pTPC < 0.6) ||
 435           (pTPC >= 0.6 && hasTOF && TMath::Abs(fPIDResponse->NumberOfSigmasTOF(track, AliPID::kProton)) < 3.0))) {
 436           isProton = kTRUE;
 437       }*/
 438       if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC(track, AliPID::kPion)) < 3.0 &&
 439           (pTPC <= 0.3 ||
 440            (pTPC > 0.3 && hasTOF && TMath::Abs(fPIDResponse->NumberOfSigmasTOF(track, AliPID::kPion)) < 3.0)))  {
 441           isPion = kTRUE;
 442       }
 443       if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC(track, AliPID::kElectron)) < 3.0 &&
 444           (pTPC <= 0.3 ||
 445            (pTPC > 0.3 && hasTOF && TMath::Abs(fPIDResponse->NumberOfSigmasTOF(track, AliPID::kElectron)) < 3.0)))  {
 446           isElectron = kTRUE;
 447       }
 448     }
 449     else if (!fMC && isV0el) {
 450       // Special treatment for V0 electrons -> Look for V0+TOF electrons
 451
 452       UInt_t status = track->GetStatus();
 453       Bool_t hasTOFout  = status&AliESDtrack::kTOFout;
 454       Bool_t hasTOFtime = status&AliESDtrack::kTIME;
 455       Bool_t hasTOF     = kFALSE;
 456       if (hasTOFout && hasTOFtime)
 457         hasTOF = kTRUE;
 458
 459       // No length check for V0's
 460       if (hasTOF && TMath::Abs(fPIDResponse->NumberOfSigmasTOF(track, AliPID::kElectron)) < 3.0) {
 461         isV0plusTOFel = kTRUE;
 462       }
 463     }
 464
 465     // MC PID, SelectSpecies, P(TPC_inner), multiplicity, deltaPrimeSpecies, eta
 466     Double_t entry[6] = { (Double_t)binMC, 0., pTPC, (Double_t)multiplicity, deltaPrimeEl, eta };
 467
 468     if (fMC)  {
 469       fhPIDdataAll->Fill(entry);
 470
 471       entry[1] = 1;
 472       entry[4] = deltaPrimeKa;
 473       fhPIDdataAll->Fill(entry);
 474
 475       entry[1] = 2;
 476       entry[4] = deltaPrimePi;
 477       fhPIDdataAll->Fill(entry);
 478
 479       entry[1] = 3;
 480       entry[4] = deltaPrimePr;
 481       fhPIDdataAll->Fill(entry);
 482     }
 483     else  {
 484       for (Int_t i = 0; i < 8; i++) {
 485         if (i == 0 && isKaon)
 486           binMC = 1;
 487         else if (i == 1 && isPion)
 488           binMC = 3;
 489         else if (i == 2 && isProton)
 490           binMC = 4;
 491         else if (i == 3 && isElectron)
 492           binMC = 0;
 493         else if (i == 4 && isV0plusTOFel)
 494           binMC = 5;
 495         else if (i == 5 && isV0el)
 496           binMC = 6;
 497         else if (i == 6 && isV0pi)
 498           binMC = 7;
 499         else if (i == 7 && isV0pr)
 500           binMC = 8;
 501         else
 502           continue;
 503
 504         entry[0] = binMC;
 505
 506         entry[1] = 0;
 507         entry[4] = deltaPrimeEl;
 508         fhPIDdataAll->Fill(entry);
 509
 510         entry[1] = 1;
 511         entry[4] = deltaPrimeKa;
 512         fhPIDdataAll->Fill(entry);
 513
 514         entry[1] = 2;
 515         entry[4] = deltaPrimePi;
 516         fhPIDdataAll->Fill(entry);
 517
 518         entry[1] = 3;
 519         entry[4] = deltaPrimePr;
 520         fhPIDdataAll->Fill(entry);
 521       }
 522     }
 523   } //track loop
 524
 525   // Post output data.
 526   PostData(1, fOutputContainer);
 527
 528   // Clear the V0 PID arrays
 529   ClearV0PIDlist();
 530 }
 531
 532 //________________________________________________________________________
 533 void AliTPCPIDEtaQA::Terminate(const Option_t *)
 534 {
 535   // Called once at the end of the query
 536   /*OLD clusterQA
 537   TObjArray* output = (TObjArray*)GetOutputData(1);
 538   if (!output)
 539     return;
 540
 541   TH1F* hPhiPrimeCutEfficiency = (TH1F*)(output->FindObject("hPhiPrimeCutEfficiency"));
 542   if (!hPhiPrimeCutEfficiency)
 543     return;
 544
 545   THnSparse* hNumClustersPhiPrimePtBeforeCut = (THnSparse*)(output->FindObject("hNumClustersPhiPrimePtBeforeCut"));
 546   if (!hNumClustersPhiPrimePtBeforeCut)
 547     return;
 548
 549   THnSparse* hNumClustersPhiPrimePtAfterCut = (THnSparse*)(output->FindObject("hNumClustersPhiPrimePtAfterCut"));
 550   if (!hNumClustersPhiPrimePtAfterCut)
 551     return;
 552
 553   TH1D* hNumEntriesBefore = hNumClustersPhiPrimePtBeforeCut->Projection(0);
 554   TH1D* hNumEntriesAfter = hNumClustersPhiPrimePtAfterCut->Projection(0);
 555
 556   for (Int_t bin = 1; bin <= hNumEntriesBefore->GetXaxis()->GetNbins(); bin++)  {
 557     Double_t numEntriesBefore = hNumEntriesBefore->GetBinContent(bin);
 558
 559     if (numEntriesBefore > 0) {
 560       Double_t numEntriesAfter = hNumEntriesAfter->GetBinContent(bin);
 561       hPhiPrimeCutEfficiency->SetBinContent(bin,  numEntriesAfter / numEntriesBefore);
 562       // Errors are 100%-correlated: Accepted tracks should be a constant fraction of all tracks in given bin.
 563       // Thus: Only take statistical fluctuation from accepted tracks as error
 564       hPhiPrimeCutEfficiency->SetBinError(bin, TMath::Sqrt(numEntriesAfter) / numEntriesBefore);
 565       //                          (numEntriesAfter / numEntriesBefore) * TMath::Sqrt(1. / numEntriesAfter + 1. / numEntriesBefore));
 566     }
 567   }
 568
 569   hPhiPrimeCutEfficiency->SetDrawOption("e");
 570   */
 571 }
 572
 573
 574 //________________________________________________________________________
 575 void AliTPCPIDEtaQA::SetUpHist(THnSparse* hist, Double_t* binsPt) const
 576 {
 577   // Sets bin limits for axes which are not standard binned and the axes titles.
 578   // MC PID, SelectSpecies, P(TPC_inner), multiplicity, deltaPrimeSpecies, eta
 579   hist->SetBinEdges(2, binsPt);
 580
 581   // Set axes titles
 582   hist->GetAxis(0)->SetTitle("MC PID");
 583   hist->GetAxis(0)->SetBinLabel(1, "e");
 584   hist->GetAxis(0)->SetBinLabel(2, "K");
 585   hist->GetAxis(0)->SetBinLabel(3, "#mu");
 586   hist->GetAxis(0)->SetBinLabel(4, "#pi");
 587   hist->GetAxis(0)->SetBinLabel(5, "p");
 588   hist->GetAxis(0)->SetBinLabel(6, "V0+TOF e");
 589   hist->GetAxis(0)->SetBinLabel(7, "V0 e");
 590   hist->GetAxis(0)->SetBinLabel(8, "V0 #pi");
 591   hist->GetAxis(0)->SetBinLabel(9, "V0 p");
 592
 593   hist->GetAxis(1)->SetTitle("Select Species");
 594   hist->GetAxis(1)->SetBinLabel(1, "e");
 595   hist->GetAxis(1)->SetBinLabel(2, "K");
 596   hist->GetAxis(1)->SetBinLabel(3, "#pi");
 597   hist->GetAxis(1)->SetBinLabel(4, "p");
 598
 599   hist->GetAxis(2)->SetTitle("p_{TPC_inner} (GeV/c)");
 600
 601   hist->GetAxis(3)->SetTitle("Event multiplicity");
 602
 603   hist->GetAxis(4)->SetTitle("TPC #Delta'{species}");
 604
 605   hist->GetAxis(5)->SetTitle("#eta");
 606
 607   //hist->Sumw2();
 608 }