- Fixed compiler warnings/errors
[u/mrichter/AliRoot.git] / 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->GetNumberOfTracks();
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] = { binMC, 0, pTPC, 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 }