Update PR task: drathee
[u/mrichter/AliRoot.git] / PWGCF / EBYE / PIDFluctuation / task / AliEbyEPidRatioDCA.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
29#include "AliESDEvent.h"
30#include "AliStack.h"
31#include "AliMCEvent.h"
32#include "AliESDtrackCuts.h"
33
34#include "AliAODEvent.h"
35#include "AliAODMCParticle.h"
36
37#include "AliEbyEPidRatioDCA.h"
38
39using namespace std;
40
41ClassImp(AliEbyEPidRatioDCA)
42
43//________________________________________________________________________
44AliEbyEPidRatioDCA::AliEbyEPidRatioDCA() :
45 AliEbyEPidRatioBase("DCA", "DCA"),
46 fESDTrackCutsBkg(NULL),
47 fHnDCA(NULL) {
48 // Constructor
49
50 AliLog::SetClassDebugLevel("AliEbyEPidRatioDCA",10);
51}
52
53//________________________________________________________________________
54AliEbyEPidRatioDCA::~AliEbyEPidRatioDCA() {
55 // Destructor
56}
57
58//________________________________________________________________________
59void AliEbyEPidRatioDCA::Process() {
60 Float_t etaRange[2];
61 fESDTrackCutsBkg->GetEtaRange(etaRange[0],etaRange[1]);
62
63 Float_t ptRange[2];
64 fESDTrackCutsBkg->GetPtRange(ptRange[0],ptRange[1]);
65
66 // -- Track Loop
67 for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
68
69 AliVTrack *track = (fESD) ? static_cast<AliVTrack*>(fESD->GetTrack(idxTrack)) : static_cast<AliVTrack*>(fAOD->GetTrack(idxTrack));
70
71 // -- Check if track is accepted for basic parameters
72 if (!fHelper->IsTrackAcceptedBasicCharged(track))
73 continue;
74
75 // -- Check if accepted - ESD
76 if (fESD && !fESDTrackCutsBkg->AcceptTrack(dynamic_cast<AliESDtrack*>(track)))
77 continue;
78
79 // -- Check if accepted - AOD
80 if (fAOD){
81 AliAODTrack * trackAOD = dynamic_cast<AliAODTrack*>(track);
82
83 if (!trackAOD) {
84 AliError("Pointer to dynamic_cast<AliAODTrack*>(track) = ZERO");
85 continue;
86 }
87 if (!trackAOD->TestFilterBit(fAODtrackCutBit))
88 continue;
89
90 // -- Check if in pT and eta range (is done in ESDTrackCuts for ESDs)
91 if(!(track->Pt() > ptRange[0] && track->Pt() <= ptRange[1] && TMath::Abs(track->Eta()) <= etaRange[1]))
92 continue;
93 }
94
95 Int_t gPdgCode = 0;
96 Int_t iPid = 0;
97 Double_t pid[3];
98 if (fHelper->IsTrackAcceptedPID(track, pid, (AliPID::kPion))) { iPid = 1; gPdgCode = 211;}
99 else if (fHelper->IsTrackAcceptedPID(track, pid, (AliPID::kKaon))) { iPid = 2; gPdgCode = 321;}
100 else if (fHelper->IsTrackAcceptedPID(track, pid, (AliPID::kProton))){ iPid = 3; gPdgCode = 2212;}
101 else iPid = 0;
102
15bb9247 103 // cout << " --- DCA ---- " << iPid << " " << gPdgCode << endl;
0a28d543 104 Double_t yP;
105 if (!fHelper->IsTrackAcceptedRapidity(track, yP, iPid))
106 continue;
107
108 Bool_t isDCArAccepted = fHelper->IsTrackAcceptedDCA(track);
109
110 // -- Check if accepted with thighter DCA cuts
111 // ?!?!? How to mimic this in AODs?
112 if (fESD && !fESDTrackCuts->AcceptTrack(dynamic_cast<AliESDtrack*>(track)))
113 isDCArAccepted = kFALSE;
114
115 // -- Check for contamination
116 Int_t contIdx = (fIsMC) ? GetContIdxTrack(TMath::Abs(track->GetLabel()), track->Charge(), gPdgCode) : 1;
117
118 // -- Get DCAs (dca_r, dca_z, sigma_xy, sigma_xy_z, sigma_z)
119 Float_t dca[2], cov[3]; //
120 if (fESD)
121 (dynamic_cast<AliESDtrack*>(track))->GetImpactParameters(dca, cov);
122 else
123 dca[0] = 1.;
124
125 // -- Fill THnSparse
126
127 if(iPid != 0) {
f7ea34d2 128 Double_t hnDCA[10] = {fCentralityBin,0,
0a28d543 129 track->Eta(),
130 yP,
131 track->Phi(),
132 track->Pt(),
133 track->Charge(),
134 contIdx,
135 isDCArAccepted,
f7ea34d2 136 dca[0]};
0a28d543 137 fHnDCA->Fill(hnDCA);
138 }
f7ea34d2 139
140 Double_t hnDCA[10] = {fCentralityBin, iPid,
141 track->Eta(),
142 yP,
143 track->Phi(),
144 track->Pt(),
145 track->Charge(),
146 contIdx,
147 isDCArAccepted,
148 dca[0]};
0a28d543 149 fHnDCA->Fill(hnDCA);
150
151
152 } // for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
153
154 return;
155}
156
157//________________________________________________________________________
158void AliEbyEPidRatioDCA::CreateHistograms() {
f7ea34d2 159 Int_t binHnDCA[11] = {AliEbyEPidRatioHelper::fgkfHistNBinsCent,4,
160 AliEbyEPidRatioHelper::fgkfHistNBinsEta,
161 AliEbyEPidRatioHelper::fgkfHistNBinsRap,
162 AliEbyEPidRatioHelper::fgkfHistNBinsPhi,
163 AliEbyEPidRatioHelper::fgkfHistNBinsPt,
164 AliEbyEPidRatioHelper::fgkfHistNBinsSign,
165 4, 2, 77};
166
167 Double_t minHnDCA[11] = {AliEbyEPidRatioHelper::fgkfHistRangeCent[0],-0.5,
168 AliEbyEPidRatioHelper::fgkfHistRangeEta[0],
169 AliEbyEPidRatioHelper::fgkfHistRangeRap[0],
170 AliEbyEPidRatioHelper::fgkfHistRangePhi[0],
171 AliEbyEPidRatioHelper::fgkfHistRangePt[0],
172 AliEbyEPidRatioHelper::fgkfHistRangeSign[0],
0a28d543 173 0.5, -0.5, -3.,-0.5};
174
f7ea34d2 175 Double_t maxHnDCA[11] = {AliEbyEPidRatioHelper::fgkfHistRangeCent[1],3.5,
176 AliEbyEPidRatioHelper::fgkfHistRangeEta[1],
177 AliEbyEPidRatioHelper::fgkfHistRangeRap[1],
178 AliEbyEPidRatioHelper::fgkfHistRangePhi[1],
179 AliEbyEPidRatioHelper::fgkfHistRangePt[1],
180 AliEbyEPidRatioHelper::fgkfHistRangeSign[1],
0a28d543 181 4.5, 1.5, 3., 3.5};
182
183
f7ea34d2 184 fHnDCA = new THnSparseD("hnDCA", "cent:pid:etaRec:yRec:phiRec:ptRec:sign:contPart:contSign:DCArAccepted:DCAr", 11, binHnDCA, minHnDCA, maxHnDCA);
0a28d543 185
186 fHnDCA->Sumw2();
f7ea34d2 187
0a28d543 188 fHnDCA->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
f7ea34d2 189 fHnDCA->GetAxis(1)->SetTitle("N_{ch}|N_{#pi}|N_{K}|N_{p}"); // 0 | 1 | 2 | 3
190 fHnDCA->GetAxis(2)->SetTitle("#eta_{Rec}"); // eta [-0.9, 0.9]
191 fHnDCA->GetAxis(3)->SetTitle("#it{y}_{Rec}"); // rapidity [-0.5, 0.5]
192 fHnDCA->GetAxis(4)->SetTitle("#varphi_{Rec} (rad)"); // phi [ 0. , 2Pi]
193 fHnDCA->GetAxis(5)->SetTitle("#it{p}_{T,Rec} (GeV/#it{c})"); // pT [ 0.2, 2.6]
194 fHnDCA->GetAxis(6)->SetTitle("sign"); // -1 | 0 | +1
195 fHnDCA->GetAxis(7)->SetTitle("contPart"); // 1 primary | 2 missId | 3 from WeakDecay | 4 p from Material
196 fHnDCA->GetAxis(8)->SetTitle("DCArAccepted"); // 0 not accepted | 1 accepted
197 fHnDCA->GetAxis(9)->SetTitle("DCAr"); // DCAr [-3, 3]
0a28d543 198
0a28d543 199
6ce4ad92 200 fHelper->BinLogAxis(fHnDCA, 5, fESDTrackCuts);
201 fHelper->BinLogAxis(fHnDCA, 5, fESDTrackCutsBkg);
0a28d543 202
203 // -- Set binning for DCAr
204 Double_t binsDCAr[77] = {-3.,-2.85,-2.7,-2.55,-2.4,-2.25,-2.1,-1.95,-1.8,-1.65,-1.5,-1.35,-1.2,-1.05,-0.9,-0.75,-0.6,-0.45,-0.3,-0.285,-0.27,-0.255,-0.24,-0.225,-0.21,-0.195,-0.18,-0.165,-0.15,-0.135,-0.12,-0.105,-0.09,-0.075,-0.06,-0.045,-0.03,-0.015,0.,0.015,0.03,0.045,0.06,0.075,0.09,0.105,0.12,0.135,0.15,0.165,0.18,0.195,0.21,0.225,0.24,0.255,0.27,0.285,0.3,0.45,0.6,0.75,0.9,1.05,1.2,1.35,1.5,1.65,1.8,1.95,2.1,2.25,2.4,2.55,2.7,2.85,3.};
f7ea34d2 205 fHnDCA->GetAxis(9)->Set(76, binsDCAr);
0a28d543 206
207 // ------------------------------------------------------------------
208
209 return;
210}
211
212//________________________________________________________________________
213Int_t AliEbyEPidRatioDCA::GetContIdxTrack(Int_t label, Int_t sign, Int_t gPdgCode) {
214 Int_t contIdx = -1;
215
216 AliVParticle* particle = (fESD) ? fMCEvent->GetTrack(label) : static_cast<AliVParticle*>(fArrayMC->At(label));
217 if (!particle)
218 return contIdx;
219
220 Bool_t isPhysicalPrimary = (fESD) ? fStack->IsPhysicalPrimary(label): (static_cast<AliAODMCParticle*>(particle))->IsPhysicalPrimary();
221 Bool_t isSecondaryFromWeakDecay = (fESD) ? fStack->IsSecondaryFromWeakDecay(label) : (static_cast<AliAODMCParticle*>(particle))->IsSecondaryFromWeakDecay();
222 Bool_t isSecondaryFromMaterial = (fESD) ? fStack->IsSecondaryFromMaterial(label) : (static_cast<AliAODMCParticle*>(particle))->IsSecondaryFromMaterial();
223
224 if (isPhysicalPrimary) {
225 if (gPdgCode == 0) {
226 // -- Check if correctly identified
227 if (particle->PdgCode() == (sign*gPdgCode))
228 contIdx = 1;
229 // -- MissIdentification
230 else
231 contIdx = 2;
232 }
233 else
234 contIdx = 1;
235 }
236
237 // -- Check if secondaries from material or weak decay
238 else if(isSecondaryFromWeakDecay)
239 contIdx = 3;
240 else if (isSecondaryFromMaterial)
241 contIdx = 4;
242 else
243 contIdx = -1;
244
245 return contIdx;
246}
247
248
249