Update PR task: drathee
[u/mrichter/AliRoot.git] / PWGCF / EBYE / PIDFluctuation / task / AliEbyEPidRatioQA.cxx
CommitLineData
0a28d543 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: ALICE Offline. *
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// AliEbyE Analysis for Particle Ratio Fluctuation //
19// Deepika Rathee | Satyajit Jena //
20// drathee@cern.ch | sjena@cern.ch //
21// Date: Wed Jul 9 18:38:30 CEST 2014 //
22// New approch to find particle ratio to reduce memory //
23// (Test Only) //
24//=========================================================================//
25
26#include "TMath.h"
27#include "TAxis.h"
28#include "TSystem.h"
29#include "TProfile.h"
30#include "TH2F.h"
31#include "TH3F.h"
32#include "TFile.h"
33#include "TPRegexp.h"
34#include "AliStack.h"
35#include "AliMCEvent.h"
36#include "AliMCParticle.h"
37#include "AliESDtrackCuts.h"
38#include "AliESDInputHandler.h"
39#include "AliESDpid.h"
40#include "AliCentrality.h"
41#include "AliTracker.h"
42#include "AliAODInputHandler.h"
43#include "AliAODEvent.h"
44#include "AliAODTrack.h"
45#include "AliAODMCParticle.h"
46#include "AliEbyEPidRatioQA.h"
47
48using namespace std;
49
50ClassImp(AliEbyEPidRatioQA)
51//________________________________________________________________________
52AliEbyEPidRatioQA::AliEbyEPidRatioQA() :
53 AliEbyEPidRatioBase("QA", "QA"),
54
f7ea34d2 55 fHnQAa(NULL), fHnQAb(NULL) {
0a28d543 56 // Constructor
57
58 AliLog::SetClassDebugLevel("AliEbyEPidRatioQA",10);
59}
60
61//________________________________________________________________________
62AliEbyEPidRatioQA::~AliEbyEPidRatioQA() {
63 // Destructor
64
65 return;
66}
67
68//________________________________________________________________________
69void AliEbyEPidRatioQA::CreateHistograms() {
70
f7ea34d2 71 Int_t binHnQAa[10] = {AliEbyEPidRatioHelper::fgkfHistNBinsCent, 4,
72 AliEbyEPidRatioHelper::fgkfHistNBinsSign,
73 AliEbyEPidRatioHelper::fgkfHistNBinsPt,
74 500, 50, 50, 50, 3};
75
76 Double_t minHnQAa[10] = {AliEbyEPidRatioHelper::fgkfHistRangeCent[0],-0.5,
77 AliEbyEPidRatioHelper::fgkfHistRangeSign[0],
78 AliEbyEPidRatioHelper::fgkfHistRangePt[0],
79 30, -10,-10, -10, -0.5};
80
81 Double_t maxHnQAa[10] = {AliEbyEPidRatioHelper::fgkfHistRangeCent[1], 3.5,
82 AliEbyEPidRatioHelper::fgkfHistRangeSign[1],
83 AliEbyEPidRatioHelper::fgkfHistRangePt[1],
84 500, 10., 10., 10., 2.5};
85
86 Int_t binHnQAb[9] = {AliEbyEPidRatioHelper::fgkfHistNBinsCent, 4,
0cb52888 87 AliEbyEPidRatioHelper::fgkfHistNBinsSign,
88 AliEbyEPidRatioHelper::fgkfHistNBinsPt,
f7ea34d2 89 AliEbyEPidRatioHelper::fgkfHistNBinsEta,
90 AliEbyEPidRatioHelper::fgkfHistNBinsRap,
91 AliEbyEPidRatioHelper::fgkfHistNBinsPhi,
f7ea34d2 92 50, 50};
0a28d543 93
f7ea34d2 94 Double_t minHnQAb[9] = {AliEbyEPidRatioHelper::fgkfHistRangeCent[0],-0.5,
95 AliEbyEPidRatioHelper::fgkfHistRangeSign[0],
0cb52888 96 AliEbyEPidRatioHelper::fgkfHistRangePt[0],
97 AliEbyEPidRatioHelper::fgkfHistRangeEta[0],
98 AliEbyEPidRatioHelper::fgkfHistRangeRap[0],
99 AliEbyEPidRatioHelper::fgkfHistRangePhi[0],
f7ea34d2 100 -10,-10};
0a28d543 101
f7ea34d2 102 Double_t maxHnQAb[9] = {AliEbyEPidRatioHelper::fgkfHistRangeCent[1], 3.5,
0cb52888 103 AliEbyEPidRatioHelper::fgkfHistRangeSign[1],
104 AliEbyEPidRatioHelper::fgkfHistRangePt[1],
105 AliEbyEPidRatioHelper::fgkfHistRangeEta[1],
106 AliEbyEPidRatioHelper::fgkfHistRangeRap[1],
107 AliEbyEPidRatioHelper::fgkfHistRangePhi[1],
108 10., 10.};
0a28d543 109
f7ea34d2 110 fHnQAa = new THnSparseF("hnQAPid", "cent:pid:sign:pt:pInner:TPCsignal:nSigmaITS:nSigmaTPC:nSigmaTOF:MCisProbe", 10, binHnQAa, minHnQAa, maxHnQAa);
111 fHnQAb = new THnSparseF("hnQADca", "cent:pid:sign:eta:y:phi:pt:DCAr:DCAz", 9, binHnQAb, minHnQAb, maxHnQAb);
112 fHnQAa->Sumw2();
113 fHnQAb->Sumw2();
0a28d543 114
f7ea34d2 115 fHnQAa->GetAxis(0)->SetTitle("centrality"); // 0-5|5-10|10-20|20-30|30-40|40-50|50-60|60-70|70-80|80-90 --> 10 bins
116 fHnQAa->GetAxis(1)->SetTitle("N_{ch}|N_{#pi}|N_{K}|N_{p}"); // 0 | 1 | 2 | 3
117 fHnQAa->GetAxis(2)->SetTitle("sign"); // -1 | 0 | +1
118 fHnQAa->GetAxis(3)->SetTitle("#it{p}_{T} (GeV/#it{c})"); // pt [ 0.46, 2.22]
119 fHnQAa->GetAxis(4)->SetTitle("#it{p}_{Inner} (GeV/#it{c})"); // pInner [ 0.1, 3.0]
120 fHnQAa->GetAxis(5)->SetTitle("TPC signal"); // TPCsignal [30, 500]
121 fHnQAa->GetAxis(6)->SetTitle("n #sigma ITS"); // nSigma ITS [-10, 10]
122 fHnQAa->GetAxis(7)->SetTitle("n #sigma TPC"); // nSigma TPC [-10, 10]
123 fHnQAa->GetAxis(8)->SetTitle("n #sigma TOF"); // nSigma TOF [-10, 10]
124 fHnQAa->GetAxis(9)->SetTitle("MCisProbe"); // 0 | 1 (isProbeParticle) | 2 (isProbeParticle wrong sign)
125
126 fHnQAb->GetAxis(0)->SetTitle("centrality"); // 0-5|5-10|10-20|20-30|30-40|40-50|50-60|60-70|70-80|80-90 --> 10 bins
127 fHnQAb->GetAxis(1)->SetTitle("N_{ch}|N_{#pi}|N_{K}|N_{p}"); // 0 | 1 | 2 | 3
128 fHnQAb->GetAxis(2)->SetTitle("sign"); // -1 | 0 | +1
0cb52888 129 fHnQAb->GetAxis(3)->SetTitle("#it{p}_{T} (GeV/#it{c})"); // pt [ 0.46, 2.22]
130 fHnQAb->GetAxis(4)->SetTitle("#eta"); // eta [-0.9, 0.9]
131 fHnQAb->GetAxis(5)->SetTitle("#it{y}"); // rapidity [-0.5, 0.5]
132 fHnQAb->GetAxis(6)->SetTitle("#varphi"); // phi [ 0. , 2Pi]
f7ea34d2 133 fHnQAb->GetAxis(7)->SetTitle("DCAr"); // DCAr [-10, 10]
134 fHnQAb->GetAxis(8)->SetTitle("DCAz"); // DCAz [-10, 10]
135 // fHnQA->GetAxis(9)->SetTitle("n #sigma #sqrt(Cdd)/DCAr"); // nSigma DCAr [-10, 10]
136 // fHnQA->GetAxis(10)->SetTitle("n #sigma #sqrt(Czz)/DCAz"); // nSigma DCAz [-10, 10]
0a28d543 137
f7ea34d2 138 fHelper->BinLogAxis(fHnQAa, 3);
0cb52888 139 fHelper->BinLogAxis(fHnQAb, 3);
f7ea34d2 140
141
0a28d543 142
f7ea34d2 143
0a28d543 144
145 return;
146}
147
148//________________________________________________________________________
149void AliEbyEPidRatioQA::Process() {
150 Float_t etaRange[2];
151 fESDTrackCuts->GetEtaRange(etaRange[0],etaRange[1]);
152
153 Float_t ptRange[2];
154 fESDTrackCuts->GetPtRange(ptRange[0],ptRange[1]);
155
156 // -- Track Loop
157 for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
15bb9247 158
0a28d543 159 AliVTrack *track = (fESD) ? static_cast<AliVTrack*>(fESD->GetTrack(idxTrack)) : static_cast<AliVTrack*>(fAOD->GetTrack(idxTrack));
160
161 // -- Check if track is accepted for basic parameters
162 if (!fHelper->IsTrackAcceptedBasicCharged(track))
163 continue;
164
165 // -- Check if accepted - ESD
15bb9247 166 if (fESD){
167 if (!fESDTrackCuts->AcceptTrack(dynamic_cast<AliESDtrack*>(track)))
168 continue;
169 }
170
171
0a28d543 172 // -- Check if accepted - AOD
173 if (fAOD){
174 AliAODTrack * trackAOD = dynamic_cast<AliAODTrack*>(track);
175
176 if (!trackAOD) {
177 AliError("Pointer to dynamic_cast<AliAODTrack*>(track) = ZERO");
178 continue;
179 }
180 if (!trackAOD->TestFilterBit(fAODtrackCutBit))
181 continue;
182
183 // -- Check if in pT and eta range (is done in ESDTrackCuts for ESDs)
184 if(!(track->Pt() > ptRange[0] && track->Pt() <= ptRange[1] && TMath::Abs(track->Eta()) <= etaRange[1]))
185 continue;
186 }
187
188 Int_t gPdgCode = 0;
189
190 Int_t iPid = 0;
191 Double_t pid[3];
192 if (fHelper->IsTrackAcceptedPID(track, pid, (AliPID::kPion))) { iPid = 1; gPdgCode = 211;}
193 else if (fHelper->IsTrackAcceptedPID(track, pid, (AliPID::kKaon))) { iPid = 2; gPdgCode = 321;}
194 else if (fHelper->IsTrackAcceptedPID(track, pid, (AliPID::kProton))){ iPid = 3; gPdgCode = 2212;}
195 else iPid = 0;
196
15bb9247 197
198 //cout << idxTrack << " --- QA ---- " << iPid << " " << gPdgCode << endl;
199
0a28d543 200 Double_t yP;
201 if ((iPid != 0) && !fHelper->IsTrackAcceptedRapidity(track, yP, iPid))
202 continue;
203
204 if (!fHelper->IsTrackAcceptedDCA(track))
205 continue;
206
207 // -- Check if probe particle
208 Int_t isProbeParticle = 0;
209 if (fIsMC) {
210 Int_t label = TMath::Abs(track->GetLabel());
211
212 AliVParticle* particle = (fESD) ? fMCEvent->GetTrack(label) : static_cast<AliVParticle*>(fArrayMC->At(label));
213 if (particle) {
214 if (TMath::Abs(particle->PdgCode()) == gPdgCode) {
215 ++isProbeParticle;
216 if (particle->PdgCode() != (track->Charge()*gPdgCode))
217 ++isProbeParticle;
218 }
219 }
220 }
221
222 // -- Get dca r/z
223 Float_t dca[] = {0.,0.}; // dca_xy, dca_z,
224 Float_t cov[] = {0.,0.,0.}; // sigma_xy, sigma_xy_z, sigma_z
225 if (fESD)
226 (static_cast<AliESDtrack*>(track))->GetImpactParameters(dca,cov);
227
f7ea34d2 228 // Float_t dcaRoverCdd = ( TMath::Sqrt(cov[0]) != 0. ) ? dca[0]/TMath::Sqrt(cov[0]) : -9.99;
229 // Float_t dcaZoverCzz = ( TMath::Sqrt(cov[2]) != 0. ) ? dca[1]/TMath::Sqrt(cov[2]) : -9.99;
0a28d543 230
231 if (iPid == 0)
232 yP = track->Eta();
233
234
235 // cout << pid[0] << " " << pid[1] << " " << pid[2] << yP << " " << iPid << endl;
236
237 if (iPid != 0) {
0a28d543 238
f7ea34d2 239 Double_t aTracka[10] = {fCentralityBin,0,track->Charge(),track->Pt(),track->GetTPCmomentum(),track->GetTPCsignal(),pid[0],pid[1],pid[2],isProbeParticle};
0cb52888 240 Double_t aTrackb[9] = {fCentralityBin,0,track->Charge(),track->Pt(),track->Eta(),yP, track->Phi(),dca[0],dca[1]};
f7ea34d2 241 fHnQAa->Fill(aTracka);
242 fHnQAb->Fill(aTrackb);
243 }
244
245 Double_t aTracka[10] = {fCentralityBin,iPid,track->Charge(),track->Pt(),track->GetTPCmomentum(),track->GetTPCsignal(),pid[0],pid[1],pid[2],isProbeParticle};
0cb52888 246 Double_t aTrackb[9] = {fCentralityBin,iPid,track->Charge(),track->Pt(),track->Eta(),yP, track->Phi(),dca[0],dca[1]};
f7ea34d2 247 fHnQAa->Fill(aTracka);
248 fHnQAb->Fill(aTrackb);
0a28d543 249
250 } // for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
251
252 return;
253}