]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGPP/pid/AliAnalysisTaskLambdaBayes.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Added an option to perform analysis w.r.t. TPC EP
[u/mrichter/AliRoot.git] / PWGPP / pid / AliAnalysisTaskLambdaBayes.cxx
CommitLineData
a8ad4709 1#include "AliAnalysisTaskLambdaBayes.h"
2
3// ROOT includes
4#include <TMath.h>
5
6// AliRoot includes
7#include "AliInputEventHandler.h"
8#include "AliAODEvent.h"
9#include "AliAODVertex.h"
10#include "AliAODTrack.h"
11#include "AliESDtrack.h"
12#include "AliCentrality.h"
13#include "AliVHeader.h"
14#include "AliAODVZERO.h"
15#include "TFile.h"
16#include "AliOADBContainer.h"
17#include "TH2F.h"
18#include "TF1.h"
19#include "AliGenHijingEventHeader.h"
20#include "AliMCEvent.h"
21#include "AliAODMCHeader.h"
22#include "AliAODMCParticle.h"
23#include "TChain.h"
24#include "AliESDtrackCuts.h"
25#include "AliESDVertex.h"
26#include "AliEventplane.h"
27#include "AliAnalysisManager.h"
28#include "TRandom.h"
29
30
31Float_t AliAnalysisTaskLambdaBayes::fPtLambdaMin[AliAnalysisTaskLambdaBayes::nPtBin] = {0.,0.25,0.5,0.75,1,1.5,2,2.5,3,3.5,4.,5.,6.};// ptmin bin
32Float_t AliAnalysisTaskLambdaBayes::fPtLambdaMax[AliAnalysisTaskLambdaBayes::nPtBin] = {0.25,0.5,0.75,1,1.5,2,2.5,3,3.5,4.,5.,6.,10.};// ptmax bin
33
34// STL includes
35//#include <iostream>
36//using namespace std;
37
38ClassImp(AliAnalysisTaskLambdaBayes)
39//_____________________________________________________________________________
40AliAnalysisTaskLambdaBayes::AliAnalysisTaskLambdaBayes():
41 AliAnalysisTaskSE(),
42 fVtxCut(10.0), // cut on |vertex| < fVtxCut
43 fEtaCut(0.8), // cut on |eta| < fEtaCut
44 fMinPt(0.15), // cut on pt > fMinPt
45 fIsMC(kFALSE),
46 fQAsw(kFALSE),
47 fNcluster(70),
48 fFilterBit(1),
49 fList(new TList()),
50 fList2(new TList()),
51 fList3(new TList()),
52 fCentrality(-1),
53 fPsi(0),
54 fPtLambdaC(0.),
55 fPhiLambdaC(0.),
56 fEtaLambda(0.),
57 fMassV0(0.),
58 fPtKp(0.),
59 fPhiKp(0.),
60 fEtaKp(0.),
61 fPtKn(0.),
62 fPhiKn(0.),
63 fEtaKn(0.),
64 fPidKp(0),
65 fPidKn(0),
66 fTOFTPCsignal(0),
67 fCombsignal(0),
68 fCutsDaughter(NULL),
69 fPIDCombined(NULL),
70 fContPid(NULL),
71 fContPid2(NULL),
520899fb 72 fContUser(NULL),
73 fContUser2(NULL),
a8ad4709 74 fNLambda(0),
75 fNpPos(0),
76 fNpNeg(0),
77 fHmismTOF(0),
9c668341 78 fHchannelTOFdistr(0),
520899fb 79 fTypeCol(2),
8c1aafae 80 fPIDuserCut(NULL),
81 fToEP(kFALSE)
a8ad4709 82{
83 // Default constructor (should not be used)
84 fList->SetName("contLambdaBayes1");
85 fList2->SetName("contLambdaBayes2");
86 fList3->SetName("contLambdaBayes3");
87
88 fList->SetOwner(kTRUE);
89 fList2->SetOwner(kTRUE);
90 fList3->SetOwner(kTRUE);
91
92 TFile *fmism = new TFile("$ALICE_ROOT/TOF/data/TOFmismatchDistr.root");
93 fHmismTOF = (TH1F *) fmism->Get("TOFmismDistr");
94
95 TFile *fchDist = new TFile("$ALICE_ROOT/TOF/data/TOFchannelDist.root");
96 fHchannelTOFdistr = (TH1D *) fchDist->Get("hTOFchanDist");
f4c4ac9e 97
98 for(Int_t i=0;i < nCentrBin;i++){
99 fElTOF[i] = NULL;
100 fElTPC[i] = NULL;
101 fPiTOF[i] = NULL;
102 fPiTPC[i] = NULL;
103 fKaTOF[i] = NULL;
104 fKaTPC[i] = NULL;
105 fPrTOF[i] = NULL;
106 fPrTPC[i] = NULL;
107 }
108 for(Int_t i=0;i < 4;i++){
109 hMatching[i] = NULL;
110 hTracking[i] = NULL;
111 }
112
113 for(Int_t i=0;i < 1000;i++){
114 fPhiLambda[i] = 0.0;
115 fPtLambda[i] = 0.0;
116 fIPPos[i] = 0;
117 fIPNeg[i] = 0;
118 fIpP[i] = 0;
119 fIpN[i] = 0;
120 fMassLambda[i] = 0.0;
121 }
a8ad4709 122}
123
124//______________________________________________________________________________
125AliAnalysisTaskLambdaBayes::AliAnalysisTaskLambdaBayes(const char *name):
126 AliAnalysisTaskSE(name),
127 fVtxCut(10.0), // cut on |vertex| < fVtxCut
128 fEtaCut(0.8), // cut on |eta| < fEtaCut
129 fMinPt(0.15), // cut on pt > fMinPt
130 fIsMC(kFALSE),
131 fQAsw(kFALSE),
132 fNcluster(70),
133 fFilterBit(1),
134 fList(new TList()),
135 fList2(new TList()),
136 fList3(new TList()),
137 fCentrality(-1),
138 fPsi(0),
139 fPtLambdaC(0.),
140 fPhiLambdaC(0.),
141 fEtaLambda(0.),
142 fMassV0(0.),
143 fPtKp(0.),
144 fPhiKp(0.),
145 fEtaKp(0.),
146 fPtKn(0.),
147 fPhiKn(0.),
148 fEtaKn(0.),
149 fPidKp(0),
150 fPidKn(0),
151 fTOFTPCsignal(0),
152 fCombsignal(0),
153 fCutsDaughter(NULL),
154 fPIDCombined(NULL),
155 fContPid(NULL),
156 fContPid2(NULL),
520899fb 157 fContUser(NULL),
158 fContUser2(NULL),
a8ad4709 159 fNLambda(0),
160 fNpPos(0),
161 fNpNeg(0),
162 fHmismTOF(0),
9c668341 163 fHchannelTOFdistr(0),
520899fb 164 fTypeCol(2),
8c1aafae 165 fPIDuserCut(NULL),
166 fToEP(kFALSE)
a8ad4709 167{
168
169 DefineOutput(1, TList::Class());
170 DefineOutput(2, TList::Class());
171 DefineOutput(3, TList::Class());
172 DefineOutput(4, TList::Class());
173
174 // Output slot #1 writes into a TTree
175 fList->SetName("contLambdaBayes1");
176 fList2->SetName("contLambdaBayes2");
177 fList3->SetName("contLambdaBayes3");
178
179 fList->SetOwner(kTRUE);
180 fList2->SetOwner(kTRUE);
181 fList3->SetOwner(kTRUE);
182
183 TFile *fmism = new TFile("$ALICE_ROOT/TOF/data/TOFmismatchDistr.root");
184 fHmismTOF = (TH1F *) fmism->Get("TOFmismDistr");
185
186 TFile *fchDist = new TFile("$ALICE_ROOT/TOF/data/TOFchannelDist.root");
187 fHchannelTOFdistr = (TH1D *) fchDist->Get("hTOFchanDist");
f4c4ac9e 188
189 for(Int_t i=0;i < nCentrBin;i++){
190 fElTOF[i] = NULL;
191 fElTPC[i] = NULL;
192 fPiTOF[i] = NULL;
193 fPiTPC[i] = NULL;
194 fKaTOF[i] = NULL;
195 fKaTPC[i] = NULL;
196 fPrTOF[i] = NULL;
197 fPrTPC[i] = NULL;
198 }
199 for(Int_t i=0;i < 4;i++){
200 hMatching[i] = NULL;
201 hTracking[i] = NULL;
202 }
203
204 for(Int_t i=0;i < 1000;i++){
205 fPhiLambda[i] = 0.0;
206 fPtLambda[i] = 0.0;
207 fIPPos[i] = 0;
208 fIPNeg[i] = 0;
209 fIpP[i] = 0;
210 fIpN[i] = 0;
211 fMassLambda[i] = 0.0;
212 }
a8ad4709 213}
214//_____________________________________________________________________________
215AliAnalysisTaskLambdaBayes::~AliAnalysisTaskLambdaBayes()
216{
217
218}
219
220//______________________________________________________________________________
221void AliAnalysisTaskLambdaBayes::UserCreateOutputObjects()
222{
223
224 fPIDCombined=new AliPIDCombined;
225 fPIDCombined->SetDefaultTPCPriors();
226 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTOF);
227
228 Float_t invMmin = 1.115-0.005*8;
229 Float_t invMmax = 1.115+0.005*8;
230
231 const Int_t nBinPid = 14;
232
8c1aafae 233 Int_t binPid[nBinPid] = {1/*ptKs*/,1+7*(!fToEP)/*EtaPiP*/,20/*pt+*/,1/*pt-*/,5/*P+*/,1/*P-*/,2/*TOFmatch+*/,1/*TOFmatch-*/,2/*istrue*/,4/*Nsigma+*/,1/*Nsigma-*/,1+4*(fToEP)/*DeltaPhi+*/,1/*DeltaPhi-*/,1+4*(fToEP)/*Psi*/};
a8ad4709 234
8c1aafae 235 Int_t binPid2[nBinPid] = {1/*ptKs*/,1+7*(!fToEP)/*EtaPiN*/,1/*pt+*/,20/*pt-*/,1/*P+*/,5/*P-*/,1/*TOFmatch+*/,2/*TOFmatch-*/,2/*istrue*/,1/*Nsigma+*/,4/*Nsigma-*/,1/*DeltaPhi+*/,1+4*(fToEP)/*DeltaPhi-*/,1+4*(fToEP)/*Psi*/};
a8ad4709 236
237 fContPid = new AliPIDperfContainer("contPID",nBinPid,binPid);
238 fContPid->SetTitleX("M_{#Lambda}");
239 fContPid->SetTitleY("centrality (%)");
240 fContPid->SetVarName(0,"p_{T}^{#Lambda}");
241 fContPid->SetVarRange(0,0,10);
242 fContPid->SetVarName(1,"#eta^{#Lambda}");
243 fContPid->SetVarRange(1,-0.8,0.8);
244 fContPid->SetVarName(2,"p_{T}^{Kp}");
245 fContPid->SetVarRange(2,0.3,4.3);
246 fContPid->SetVarName(3,"p_{T}^{Kn}");
247 fContPid->SetVarRange(3,0.3,4.3);
248 fContPid->SetVarName(4,"BayesProb^{Kp}");
249 fContPid->SetVarRange(4,0,1.);
250 fContPid->SetVarName(5,"BayesProb^{Kn}");
251 fContPid->SetVarRange(5,0,1.);
252 fContPid->SetVarName(6,"isTOFmatch^{Kp}");
253 fContPid->SetVarRange(6,-0.5,1.5);
254 fContPid->SetVarName(7,"isTOFmatch^{Kn}");
255 fContPid->SetVarRange(7,-0.5,1.5);
256 fContPid->SetVarName(8,"isLambdaTrue^{Kn}");
257 fContPid->SetVarRange(8,-0.5,1.5);
258 fContPid->SetVarName(9,"N#sigma^{Kp}");
259 fContPid->SetVarRange(9,1.25,6.25);
260 fContPid->SetVarName(10,"N#sigma^{Kn}");
261 fContPid->SetVarRange(10,1.25,6.25);
262 fContPid->SetVarName(11,"#Delta#phi^{Kp}");
263 fContPid->SetVarRange(11,-TMath::Pi(),TMath::Pi());
264 fContPid->SetVarName(12,"#Delta#phi^{Kn}");
265 fContPid->SetVarRange(12,-TMath::Pi(),TMath::Pi());
266 fContPid->SetVarName(13,"#Psi");
267 fContPid->SetVarRange(13,-TMath::Pi()/2,TMath::Pi()/2);
268
269 fContPid2 = new AliPIDperfContainer("contPID2",nBinPid,binPid2);
270 fContPid2->SetTitleX("M_{#Lambda}");
271 fContPid2->SetTitleY("centrality (%)");
272 fContPid2->SetVarName(0,"p_{T}^{#Lambda}");
273 fContPid2->SetVarRange(0,0,10);
274 fContPid2->SetVarName(1,"#eta^{#Lambda}");
275 fContPid2->SetVarRange(1,-0.8,0.8);
276 fContPid2->SetVarName(2,"p_{T}^{Kp}");
277 fContPid2->SetVarRange(2,0.3,4.3);
278 fContPid2->SetVarName(3,"p_{T}^{Kn}");
279 fContPid2->SetVarRange(3,0.3,4.3);
280 fContPid2->SetVarName(4,"BayesProb^{Kp}");
281 fContPid2->SetVarRange(4,0,1.);
282 fContPid2->SetVarName(5,"BayesProb^{Kn}");
283 fContPid2->SetVarRange(5,0,1.);
284 fContPid2->SetVarName(6,"isTOFmatch^{Kp}");
285 fContPid2->SetVarRange(6,-0.5,1.5);
286 fContPid2->SetVarName(7,"isTOFmatch^{Kn}");
287 fContPid2->SetVarRange(7,-0.5,1.5);
288 fContPid2->SetVarName(8,"isLambdaTrue^{Kn}");
289 fContPid2->SetVarRange(8,-0.5,1.5);
290 fContPid2->SetVarName(9,"N#sigma^{Kp}");
291 fContPid2->SetVarRange(9,1.25,6.25);
292 fContPid2->SetVarName(10,"N#sigma^{Kn}");
293 fContPid2->SetVarRange(10,1.25,6.25);
294 fContPid2->SetVarName(11,"#Delta#phi^{Kp}");
295 fContPid2->SetVarRange(11,-TMath::Pi(),TMath::Pi());
296 fContPid2->SetVarName(12,"#Delta#phi^{Kn}");
297 fContPid2->SetVarRange(12,-TMath::Pi(),TMath::Pi());
298 fContPid2->SetVarName(13,"#Psi");
299 fContPid2->SetVarRange(13,-TMath::Pi()/2,TMath::Pi()/2);
300
301
302
303 const Int_t nDETsignal = 100; // mass
304 Double_t binDETsignal[nDETsignal+1];
305 for(Int_t i=0;i<nDETsignal+1;i++){
306 binDETsignal[i] = invMmin + i*(invMmax - invMmin) / nDETsignal;
307 }
308 const Int_t nDETsignal2 = 10; // centrality
309 Double_t binDETsignal2[nDETsignal2+1];
310 for(Int_t i=0;i<nDETsignal2+1;i++){
311 binDETsignal2[i] = i*100./ nDETsignal2;
312 }
313 fContPid->AddSpecies("Lambda",nDETsignal,binDETsignal,nDETsignal2,binDETsignal2);
314 fContPid2->AddSpecies("Lambda2",nDETsignal,binDETsignal,nDETsignal2,binDETsignal2);
315
316 fList->Add(fContPid);
317 fList->Add(fContPid2);
318
520899fb 319 const Int_t nBinUser = 6;
320 Int_t binUser[nBinUser] = {8/*Eta*/,20/*pt*/,2/*istrue*/,4/*whatSelection*/,1/*DeltaPhi*/,1/*Psi*/};
321 fContUser = new AliPIDperfContainer("contUserPID",nBinUser,binUser);
322 fContUser->SetTitleX("M_{#Lambda}");
323 fContUser->SetTitleY("centrality (%)");
324 fContUser->SetVarName(0,"#eta");
325 fContUser->SetVarRange(0,-0.8,0.8);
326 fContUser->SetVarName(1,"p_{T}");
327 fContUser->SetVarRange(1,0.3,4.3);
328 fContUser->SetVarName(2,"isLambdaTrue");
329 fContUser->SetVarRange(2,-0.5,1.5);
330 fContUser->SetVarName(3,"whatSelected"); // 0=no, 1=pi, 2=K, 3=p
331 fContUser->SetVarRange(3,-0.5,3.5);
332 fContUser->SetVarName(4,"#Delta#phi");
333 fContUser->SetVarRange(4,-TMath::Pi(),TMath::Pi());
334 fContUser->SetVarName(5,"#Psi");
335 fContUser->SetVarRange(5,-TMath::Pi()/2,TMath::Pi()/2);
336
337 fContUser2 = new AliPIDperfContainer("contUserPID2",nBinUser,binUser);
338 fContUser2->SetTitleX("M_{#Lambda}");
339 fContUser2->SetTitleY("centrality (%)");
340 fContUser2->SetVarName(0,"#eta");
341 fContUser2->SetVarRange(0,-0.8,0.8);
342 fContUser2->SetVarName(1,"p_{T}");
343 fContUser2->SetVarRange(1,0.3,4.3);
344 fContUser2->SetVarName(2,"isLambdaTrue");
345 fContUser2->SetVarRange(2,-0.5,1.5);
346 fContUser2->SetVarName(3,"whatSelected");
347 fContUser2->SetVarRange(3,-0.5,3.5);
348 fContUser2->SetVarName(4,"#Delta#phi");
349 fContUser2->SetVarRange(4,-TMath::Pi(),TMath::Pi());
350 fContUser2->SetVarName(5,"#Psi");
351 fContUser2->SetVarRange(5,-TMath::Pi()/2,TMath::Pi()/2);
352
353 fContUser->AddSpecies("Lambda",nDETsignal,binDETsignal,nDETsignal2,binDETsignal2);
354 fContUser2->AddSpecies("Lambda2",nDETsignal,binDETsignal,nDETsignal2,binDETsignal2);
355
356 fList->Add(fContUser);
357 fList->Add(fContUser2);
358
a8ad4709 359 hMatching[0] = new TH2F("hMatchAll","TOF matched (all);p_{T} (GeV/#it{c});centrality (%)",50,0,10,nDETsignal2,0,100);
360 hMatching[1] = new TH2F("hMatchPi","TOF matched (#pi);p_{T} (GeV/#it{c});centrality (%)",50,0,10,nDETsignal2,0,100);
361 hMatching[2] = new TH2F("hMatchKa","TOF matched (K);p_{T} (GeV/#it{c});centrality (%)",50,0,10,nDETsignal2,0,100);
362 hMatching[3] = new TH2F("hMatchPr","TOF matched (p);p_{T} (GeV/#it{c});centrality (%)",50,0,10,nDETsignal2,0,100);
363
364 hTracking[0] = new TH2F("hTrackingAll","TPC tracks (all);p_{T} (GeV/#it{c});centrality (%)",50,0,10,nDETsignal2,0,100);
365 hTracking[1] = new TH2F("hTrackingPi","TPC tracks (#pi);p_{T} (GeV/#it{c});centrality (%)",50,0,10,nDETsignal2,0,100);
366 hTracking[2] = new TH2F("hTrackingKa","TPC tracks (K);p_{T} (GeV/#it{c});centrality (%)",50,0,10,nDETsignal2,0,100);
367 hTracking[3] = new TH2F("hTrackingPr","TPC tracks (p);p_{T} (GeV/#it{c});centrality (%)",50,0,10,nDETsignal2,0,100);
368
369 fList2->Add(hMatching[0]);
370 fList2->Add(hMatching[1]);
371 fList2->Add(hMatching[2]);
372 fList2->Add(hMatching[3]);
373 fList2->Add(hTracking[0]);
374 fList2->Add(hTracking[1]);
375 fList2->Add(hTracking[2]);
376 fList2->Add(hTracking[3]);
377
378 fTOFTPCsignal = new TH2F("hTOFTPCsignal","TOF-TPC signal for protons (1.2 < p_{T} < 1.3 GeV/#it{c}, cent. = 0-20);N#sigma_{TOF};N#sigma_{TPC}",100,-5,5,100,-5,5);
379 fList2->Add(fTOFTPCsignal);
380 fCombsignal = new TH1F("hCombsignal","Combined signal for protons (1.2 < p_{T} < 1.3 GeV/#it{c}, cent. = 0-20);N#sigma",100,0,5);
381 fList2->Add(fCombsignal);
382
383 // QA plots
384 char name[100],title[400];
385 for(Int_t i=0;i < nCentrBin;i++){
386 snprintf(name,100,"hPiTPCQA_%i",i);
387 snprintf(title,400,"TPC signal for #pi cent=%i-%i%c;p_{T} (GeV/#it{c});N#sigma",i*10,(i+1)*10,'%');
388 fPiTPC[i] = new TH2F(name,title,50,0,4,200,-10,10);
389 fList3->Add(fPiTPC[i]);
390
391 snprintf(name,100,"hKaTPCQA_%i",i);
392 snprintf(title,400,"TPC signal for K cent=%i-%i%c;p_{T} (GeV/#it{c});N#sigma",i*10,(i+1)*10,'%');
393 fKaTPC[i] = new TH2F(name,title,50,0,4,200,-10,10);
394 fList3->Add(fKaTPC[i]);
395
396 snprintf(name,100,"hPrTPCQA_%i",i);
397 snprintf(title,400,"TPC signal for p cent=%i-%i%c;p_{T} (GeV/#it{c});N#sigma",i*10,(i+1)*10,'%');
398 fPrTPC[i] = new TH2F(name,title,50,0,4,200,-10,10);
399 fList3->Add(fPrTPC[i]);
400
401 snprintf(name,100,"hElTPCQA_%i",i);
402 snprintf(title,400,"TPC signal for e cent=%i-%i%c;p_{T} (GeV/#it{c});N#sigma",i*10,(i+1)*10,'%');
403 fElTPC[i] = new TH2F(name,title,50,0,4,200,-10,10);
404 fList3->Add(fElTPC[i]);
405
406 snprintf(name,100,"hPiTOFQA_%i",i);
407 snprintf(title,400,"TOF signal for #pi cent=%i-%i%c;p_{T} (GeV/#it{c});N#sigma",i*10,(i+1)*10,'%');
408 fPiTOF[i] = new TH2F(name,title,50,0,4,200,-10,10);
409 fList3->Add(fPiTOF[i]);
410
411 snprintf(name,100,"hKaTOFQA_%i",i);
412 snprintf(title,400,"TOF signal for K cent=%i-%i%c;p_{T} (GeV/#it{c});N#sigma",i*10,(i+1)*10,'%');
413 fKaTOF[i] = new TH2F(name,title,50,0,4,200,-10,10);
414 fList3->Add(fKaTOF[i]);
415
416 snprintf(name,100,"hPrTOFQA_%i",i);
417 snprintf(title,400,"TOF signal for p cent=%i-%i%c;p_{T} (GeV/#it{c});N#sigma",i*10,(i+1)*10,'%');
418 fPrTOF[i] = new TH2F(name,title,50,0,4,200,-10,10);
419 fList3->Add(fPrTOF[i]);
420
421 snprintf(name,100,"hElTOFQA_%i",i);
422 snprintf(title,400,"TOF signal for e cent=%i-%i%c;p_{T} (GeV/#it{c});N#sigma",i*10,(i+1)*10,'%');
423 fElTOF[i] = new TH2F(name,title,50,0,4,200,-10,10);
424 fList3->Add(fElTOF[i]);
425
426 }
427
428 // Post output data.
429 PostData(1, fList);
430 PostData(2, fList2);
431 PostData(3, fList3);
432}
433
434//______________________________________________________________________________
435void AliAnalysisTaskLambdaBayes::UserExec(Option_t *)
436{
437 // Main loop
438 // Called for each event
439
440 fOutputAOD = dynamic_cast<AliAODEvent*>(InputEvent());
441 if(!fOutputAOD){
442 Printf("%s:%d AODEvent not found in Input Manager",(char*)__FILE__,__LINE__);
443 this->Dump();
444 return;
445 }
446
447 Float_t zvtx = GetVertex(fOutputAOD);
448
449
450
451 //Get the MC object
452 if(fIsMC){
453 AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(fOutputAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
454 if (!mcHeader) {
455 AliError("Could not find MC Header in AOD");
456 return;
457 }
458 }
459
460 if (TMath::Abs(zvtx) < fVtxCut) {
461
462 SelectLambda();
463
464 //Centrality
465 Float_t v0Centr = -10.;
466 Float_t trkCentr = -10.;
467 AliCentrality *centrality = fOutputAOD->GetCentrality();
468 if (centrality){
469 trkCentr = centrality->GetCentralityPercentile("V0M");
470 v0Centr = centrality->GetCentralityPercentile("TRK");
471 }
472
cd71f9f4 473 if(!fTypeCol){
474 v0Centr=100./(fOutputAOD->GetNumberOfTracks()/12.+1);
475 trkCentr=v0Centr;
476 }
477
9c668341 478 if((TMath::Abs(v0Centr - trkCentr) < 5.0 || (fTypeCol!=2)) && v0Centr>0){ // consistency cut on centrality selection
a8ad4709 479 fCentrality = v0Centr;
480 Analyze(fOutputAOD); // Do analysis!!!!
481
482 }
483 }
484
485}
486
487//________________________________________________________________________
488void AliAnalysisTaskLambdaBayes::Analyze(AliAODEvent* aodEvent)
489{
490
491 Int_t ntrack = aodEvent->GetNumberOfTracks();
492
493 fPtKp=0.,fPhiKp=0.,fEtaKp=0.;
494 fPtKn=0.,fPhiKn=0.,fEtaKn=0.;
495 fPidKp=0,fPidKn=0;
496 fMassV0=-1;
497
498 TLorentzVector ksV;
499
500 Double_t px,py,pz,E;
501
502 Float_t invMmin = 1.115-0.005*8;
503 Float_t invMmax = 1.115+0.005*8;
504
505 Int_t icentr = 8;
506 if(fCentrality < 0) icentr = 8;
507 else if(fCentrality < 10) icentr = 0;
508 else if(fCentrality < 20) icentr = 1;
509 else if(fCentrality < 30) icentr = 2;
510 else if(fCentrality < 40) icentr = 3;
511 else if(fCentrality < 50) icentr = 4;
512 else if(fCentrality < 60) icentr = 5;
513 else if(fCentrality < 70) icentr = 6;
514 else if(fCentrality < 80) icentr = 7;
515
516 Float_t addMismatchForMC = 0.005;
faa6d35f 517 if(fCentrality < 50) addMismatchForMC += 0.005;
a8ad4709 518 if(fCentrality < 20) addMismatchForMC += 0.02;
519
faa6d35f 520 if(fTypeCol == 0 || fTypeCol == 1) addMismatchForMC = 0.005;
9c668341 521
a8ad4709 522 fPsi = 0;
523 /* Compute TPC EP */
524 Double_t Qx2 = 0, Qy2 = 0;
525 Double_t Qx3 = 0, Qy3 = 0;
526 for(Int_t iT = 0; iT < ntrack; iT++) {
527 AliAODTrack* aodTrack = aodEvent->GetTrack(iT);
528
529 if (!aodTrack){
530 continue;
531 }
532
533 Bool_t trkFlag = aodTrack->TestFilterBit(fFilterBit);
534
535 if ((TMath::Abs(aodTrack->Eta()) > 0.8) || (aodTrack->Pt() < 0.2) || (aodTrack->GetTPCNcls() < fNcluster) || !trkFlag)
536 continue;
537
538 Double_t b[2] = {-99., -99.};
539 Double_t bCov[3] = {-99., -99., -99.};
540 if (!aodTrack->PropagateToDCA(aodEvent->GetPrimaryVertex(), aodEvent->GetMagneticField(), 100., b, bCov))
541 continue;
542
543 if ((TMath::Abs(b[0]) > 3.0) || (TMath::Abs(b[1]) > 2.4))
544 continue;
545
546 Qx2 += TMath::Cos(2*aodTrack->Phi());
547 Qy2 += TMath::Sin(2*aodTrack->Phi());
548 Qx3 += TMath::Cos(3*aodTrack->Phi());
549 Qy3 += TMath::Sin(3*aodTrack->Phi());
550
551 }
552
553 fPsi = TMath::ATan2(Qy2, Qx2)/2.;
554
555 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
556 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
557 AliPIDResponse *PIDResponse=inputHandler->GetPIDResponse();
558
cd71f9f4 559// PIDResponse->GetTOFResponse().SetTrackParameter(0,0.);
560// PIDResponse->GetTOFResponse().SetTrackParameter(1,0.);
561// PIDResponse->GetTOFResponse().SetTrackParameter(2,0.018);
562// PIDResponse->GetTOFResponse().SetTrackParameter(3,50.0);
a8ad4709 563
564 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC|AliPIDResponse::kDetTOF);
565
f4c4ac9e 566
567 Double_t probP[10] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
568 Double_t probN[10] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
a8ad4709 569 Double_t nSigmaTPC,nSigmaTOF=6,nSigmaTPC2,nSigmaTOF2=6,nSigmaComb,nSigmaComb2;
570
571
572 AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
573 TClonesArray *mcArray = NULL;
574 if (mcHeader)
575 mcArray = (TClonesArray*)aodEvent->GetList()->FindObject(AliAODMCParticle::StdBranchName());
576
577 Int_t nmc = 0;
578 if(mcArray)
579 nmc = mcArray->GetEntries();
580
581 for(Int_t i=0;i < ntrack;i++){ // loop on tracks
582 AliAODTrack* track = aodEvent->GetTrack(i);
583
584 AliAODMCParticle *mcp = NULL;
585 Int_t pdg = 0;
586
587 if (!track){
588 continue;
589 }
590
591 Int_t tofMatch = (track->GetStatus() & AliVTrack::kTOFout) && (track->GetStatus() & AliVTrack::kTIME);
592
593 Int_t label = -1;
594 if(mcArray){
595 label = track->GetLabel();
596 if(label != -1 && label < nmc){
597 label = TMath::Abs(label);
598 mcp = (AliAODMCParticle*)mcArray->At(label);
599 pdg = TMath::Abs(mcp->GetPdgCode());
600 }
601 else
602 label = -1;
603 }
604 else{
605 /*UInt_t detUsed =*/ fPIDCombined->ComputeProbabilities(track, PIDResponse, probP);
606 }
607
608 if(track->TestFilterBit(fFilterBit) && TMath::Abs(track->Eta()) < 0.8 && track->Charge() > 0){
609 hTracking[0]->Fill(track->Pt(),fCentrality);
610 if(pdg == 211)
611 hTracking[1]->Fill(track->Pt(),fCentrality);
612 else if(pdg == 321)
613 hTracking[2]->Fill(track->Pt(),fCentrality);
614 else if(pdg == 2212)
615 hTracking[3]->Fill(track->Pt(),fCentrality);
616 else if(! mcp){ // Fill matching histo with the prob
617 hTracking[1]->Fill(track->Pt(),fCentrality,probP[2]);
618 hTracking[2]->Fill(track->Pt(),fCentrality,probP[3]);
619 hTracking[3]->Fill(track->Pt(),fCentrality,probP[4]);
620 }
621 }
622
623 if(!tofMatch) continue;
624
625 if(track->TestFilterBit(fFilterBit) && TMath::Abs(track->Eta()) < 0.8 && track->Charge() > 0){
626 hMatching[0]->Fill(track->Pt(),fCentrality);
627 if(pdg == 211)
628 hMatching[1]->Fill(track->Pt(),fCentrality);
629 else if(pdg == 321)
630 hMatching[2]->Fill(track->Pt(),fCentrality);
631 else if(pdg == 2212)
632 hMatching[3]->Fill(track->Pt(),fCentrality);
633 else if(! mcp){ // Fill matching histo with the prob
634 hMatching[1]->Fill(track->Pt(),fCentrality,probP[2]);
635 hMatching[2]->Fill(track->Pt(),fCentrality,probP[3]);
636 hMatching[3]->Fill(track->Pt(),fCentrality,probP[4]);
637 }
638 }
639 }
640
007df888 641// Int_t pdg1 = -1;
642// Int_t pdg2 = -1;
a8ad4709 643
644
645 // start analysis Lambda
646 for(Int_t i=0;i < ntrack;i++){ // loop on proton candidate tracks
647 AliAODTrack* KpTrack = aodEvent->GetTrack(i);
648
649 if (!KpTrack){
650 continue;
651 }
652
653 if(!(KpTrack->Pt() > 0.3 && TMath::Abs(KpTrack->Eta()) < 0.8)) continue;
654
655 Bool_t isLambda = (KpTrack->Charge() > 0);
656
657 nSigmaComb=5;
e34b28fe 658 nSigmaTOF=5;
a8ad4709 659 fPtKp=KpTrack->Pt(),fPhiKp=KpTrack->Phi(),fEtaKp=KpTrack->Eta();
660 fPidKp=0;
661
662 UInt_t detUsedP = fPIDCombined->ComputeProbabilities(KpTrack, PIDResponse, probP);
663
664 Double_t oldpP[10];
665 fPIDCombined->GetPriors(KpTrack, oldpP, PIDResponse, detUsedP);
666
667
a8ad4709 668 nSigmaTPC = PIDResponse->NumberOfSigmasTPC(KpTrack,AliPID::kProton);
a8ad4709 669
670 if(! (TMath::Abs(nSigmaTPC) < 5)) continue;
671
672 Int_t tofMatch1 = (KpTrack->GetStatus() & AliVTrack::kTOFout) && (KpTrack->GetStatus() & AliVTrack::kTIME);
9c668341 673 /*
a8ad4709 674 if(mcArray){
675 Int_t labelK = TMath::Abs(KpTrack->GetLabel());
676 AliAODMCParticle *mcp1 = (AliAODMCParticle*)mcArray->At(labelK);
9c668341 677 pdg1 = TMath::Abs(mcp1->GetPdgCode());
a8ad4709 678 }
9c668341 679 */
a8ad4709 680
681 fPidKp = Int_t(probP[4]*100);
682
683 if(tofMatch1){
684 if(!IsChannelValid(TMath::Abs(KpTrack->Eta()))){
685 // remove this tof hit
686 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC);
687 detUsedP = fPIDCombined->ComputeProbabilities(KpTrack, PIDResponse, probP);
688 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC|AliPIDResponse::kDetTOF);
689 fPidKp = Int_t(probP[4]*100);
690 tofMatch1=0;
691 }
692 else{
693 if(probP[4] > probP[3] && probP[4] > probP[2] && probP[4] > probP[0]) fPidKp += 128; // max prob
694
cd71f9f4 695 if(isLambda){
696 nSigmaTOF = PIDResponse->NumberOfSigmasTOF(KpTrack,AliPID::kElectron);
697 if(TMath::Abs(PIDResponse->NumberOfSigmasTPC(KpTrack,AliPID::kElectron))<1) fElTOF[icentr]->Fill(fPtKp,nSigmaTOF);
698 if(TMath::Abs(nSigmaTOF)<1) fElTPC[icentr]->Fill(fPtKp,PIDResponse->NumberOfSigmasTPC(KpTrack,AliPID::kElectron));
699 nSigmaTOF = PIDResponse->NumberOfSigmasTOF(KpTrack,AliPID::kPion);
700 if(TMath::Abs(PIDResponse->NumberOfSigmasTPC(KpTrack,AliPID::kPion))<1) fPiTOF[icentr]->Fill(fPtKp,nSigmaTOF);
701 if(TMath::Abs(nSigmaTOF)<1) fPiTPC[icentr]->Fill(fPtKp,PIDResponse->NumberOfSigmasTPC(KpTrack,AliPID::kPion));
702 nSigmaTOF = PIDResponse->NumberOfSigmasTOF(KpTrack,AliPID::kKaon);
703 if(TMath::Abs(PIDResponse->NumberOfSigmasTPC(KpTrack,AliPID::kKaon))<1) fKaTOF[icentr]->Fill(fPtKp,nSigmaTOF);
704 if(TMath::Abs(nSigmaTOF)<1) fKaTPC[icentr]->Fill(fPtKp,PIDResponse->NumberOfSigmasTPC(KpTrack,AliPID::kKaon));
705 }
a8ad4709 706 nSigmaTOF = PIDResponse->NumberOfSigmasTOF(KpTrack,AliPID::kProton);
cd71f9f4 707 if(isLambda){
708 if(TMath::Abs(PIDResponse->NumberOfSigmasTPC(KpTrack,AliPID::kProton))<1) fPrTOF[icentr]->Fill(fPtKp,nSigmaTOF);
709 if(TMath::Abs(nSigmaTOF)<1) fPrTPC[icentr]->Fill(fPtKp,PIDResponse->NumberOfSigmasTPC(KpTrack,AliPID::kProton));
710 }
a8ad4709 711
a8ad4709 712 if(fIsMC){
713 Float_t mismAdd = addMismatchForMC;
714 if(KpTrack->Pt() < 1) mismAdd = addMismatchForMC/KpTrack->Pt();
715
716 if(gRandom->Rndm() < mismAdd){
717 Float_t etaAbs = TMath::Abs(KpTrack->Eta());
718 Int_t channel = Int_t(4334.09 - 4758.36 * etaAbs -1989.71 * etaAbs*etaAbs + 1957.62*etaAbs*etaAbs*etaAbs);
719 channel = channel % 8736;
720 Float_t distIP = fHchannelTOFdistr->GetBinContent(channel);
721
722 // generate random time
a16898fb 723 Float_t timeRandom = fHmismTOF->GetRandom() + distIP*3.35655419905265973e+01;
a8ad4709 724 Double_t times[10];
725 KpTrack->GetIntegratedTimes(times);
726 nSigmaTOF = TMath::Abs(timeRandom - times[4])/PIDResponse->GetTOFResponse().GetExpectedSigma(KpTrack->P(), times[3], AliPID::kProton);
727 }
728 }
729
730 if(fCentrality < 20 && KpTrack->Pt() < 1.3 && KpTrack->Pt() > 1.2)fTOFTPCsignal->Fill(nSigmaTOF,nSigmaTPC);
63d1d548 731 nSigmaTOF = TMath::Abs(nSigmaTOF);
a8ad4709 732
733 if(nSigmaTOF < 2) fPidKp += 256;
734 else if(nSigmaTOF < 3) fPidKp += 512;
735 }
736 }
737
738 if(tofMatch1){
739 nSigmaComb = TMath::Sqrt(nSigmaTOF*nSigmaTOF + nSigmaTPC*nSigmaTPC);
740 if(nSigmaTOF < 5 && fCentrality < 20 && KpTrack->Pt() < 1.3 && KpTrack->Pt() > 1.2){
741 fCombsignal->Fill(nSigmaComb);
742 }
743 } else {
744 nSigmaComb = TMath::Abs(nSigmaTPC);
745 }
746
747 // use sigmaTOF instead of sigmaComb
e34b28fe 748 if(tofMatch1) nSigmaComb = nSigmaTOF;
a8ad4709 749
750 if(nSigmaComb < 2) nSigmaComb = 2;
751 else if(nSigmaComb < 3) nSigmaComb = 3;
752 else if(nSigmaComb < 5) nSigmaComb = 4.99;
753 else nSigmaComb = 6;
754
755 Int_t iks=-1;
756 for(Int_t k=0;k < fNLambda;k++){ // find the Lambda which contains the positive track
757 if(i == fIpP[k]) iks = k;
758 }
759
760 if(fPtKp > 4.299) fPtKp = 4.299;
761
762 if(iks > -1 && fIpN[iks] > -1){
763 //for(Int_t j=0;j < ntrack;j++){ // loop on negative tracks
764 Int_t j = fIpN[iks];
765 AliAODTrack* KnTrack = aodEvent->GetTrack(j);
766
767 if (!KnTrack){
768 continue;
769 }
770
e34b28fe 771 nSigmaComb2 = 5;
772 nSigmaTOF2 = 5;
773
a8ad4709 774 if(!(KnTrack->Pt() > 0.3 && TMath::Abs(KnTrack->Eta()) < 0.8)) continue;
775
776 fPtKn=KnTrack->Pt(),fPhiKn=KnTrack->Phi(),fEtaKn=KnTrack->Eta();
777 fPidKn=0;
778
779 UInt_t detUsedN = fPIDCombined->ComputeProbabilities(KnTrack, PIDResponse, probN);
780 Double_t oldpN[10];
781 fPIDCombined->GetPriors(KnTrack, oldpN, PIDResponse, detUsedN);
782
783 nSigmaTPC2 = PIDResponse->NumberOfSigmasTPC(KnTrack,AliPID::kPion);
784
785 if(! (TMath::Abs(nSigmaTPC2) < 5)) continue;
786
a8ad4709 787 Int_t tofMatch2 = (KnTrack->GetStatus() & AliVTrack::kTOFout) && (KnTrack->GetStatus() & AliVTrack::kTIME);
9c668341 788 /*
a8ad4709 789 if(mcArray){
790 Int_t labelK = TMath::Abs(KnTrack->GetLabel());
791 AliAODMCParticle *mcp2 = (AliAODMCParticle*)mcArray->At(labelK);
9c668341 792 pdg2 = TMath::Abs(mcp2->GetPdgCode());
793 }
794 */
a8ad4709 795
796 fPidKn = Int_t(probN[2]*100);
797
798 if(tofMatch2){
799 if(probN[2] > probN[3] && probN[2] > probN[4] && probN[2] > probN[0]) fPidKn += 128; // max prob
800
801 nSigmaTOF2 = PIDResponse->NumberOfSigmasTOF(KnTrack,AliPID::kPion);
802
803 nSigmaTOF2 = TMath::Abs(nSigmaTOF2);
804
805 if(nSigmaTOF2 < 2) fPidKn += 256;
806 else if(nSigmaTOF2 < 3) fPidKn += 512;
807 }
808
809 px = KpTrack->Px() + KnTrack->Px();
810 py = KpTrack->Py() + KnTrack->Py();
811 pz = KpTrack->Pz() + KnTrack->Pz();
812 E = TMath::Sqrt(KpTrack->P()*KpTrack->P() + 0.938e-01*0.938e-01);
813 E += TMath::Sqrt(KnTrack->P()*KnTrack->P()+ 1.39e-01*1.39e-01);
814
815 ksV.SetPxPyPzE(px,py,pz,E);
816 fMassV0 = fMassLambda[iks];
817
818 if(fMassV0 < invMmin || fMassV0 > invMmax) continue;
819
820
821 fPtLambdaC = ksV.Pt();
822 fEtaLambda = ksV.Eta();
823 fPhiLambdaC = ksV.Phi();
824
825 if(tofMatch2){
826 nSigmaComb2 = TMath::Sqrt(nSigmaTOF2*nSigmaTOF2+ nSigmaTPC2*nSigmaTPC2);
827 } else {
828 nSigmaComb2 = TMath::Abs(nSigmaTPC2);
829 }
830
831 // use sigmaTOF instead of sigmaComb
e34b28fe 832 if(tofMatch2) nSigmaComb2 = nSigmaTOF2;
a8ad4709 833
834 if(nSigmaComb2 < 2) nSigmaComb2 = 2;
835 else if(nSigmaComb2 < 3) nSigmaComb2 = 3;
836 else if(nSigmaComb2 < 5) nSigmaComb2 = 4.99;
837 else nSigmaComb2 = 6;
838
839 Bool_t isTrue = kFALSE;
840
841 if(mcArray){
842 Int_t labelP = TMath::Abs(KpTrack->GetLabel());
843 Int_t labelN = TMath::Abs(KnTrack->GetLabel());
844
845 if(labelP > -1 && labelN > -1){
846 AliAODMCParticle *partP = (AliAODMCParticle*)mcArray->At(labelP);
847 AliAODMCParticle *partN = (AliAODMCParticle*)mcArray->At(labelN);
848
849 Int_t mP = partP->GetMother();
850 Int_t mN = partN->GetMother();
851 if(mP == mN && mP > -1){
852 AliAODMCParticle *partM = (AliAODMCParticle*)mcArray->At(mP);
853 Int_t pdgM = partM->GetPdgCode();
854 if(TMath::Abs(pdgM) == 3122) isTrue = kTRUE;
855 }
856 }
857
858 }
859
860 Double_t deltaphi1 = KpTrack->Phi() - fPsi;
861 Double_t deltaphi2 = KnTrack->Phi() - fPsi;
862
863 while(deltaphi1 > TMath::Pi()) deltaphi1 -= TMath::Pi()*2;
864 while(deltaphi1 < -TMath::Pi()) deltaphi1 += TMath::Pi()*2;
865 while(deltaphi2 > TMath::Pi()) deltaphi2 -= TMath::Pi()*2;
866 while(deltaphi2 < -TMath::Pi()) deltaphi2 += TMath::Pi()*2;
867
868 if(fPtKn > 4.299) fPtKn = 4.299;
869
870 Float_t xTOfill[] = {fPtLambdaC,KpTrack->Eta(),fPtKp,fPtKn,(fPidKp%128)*0.01,(fPidKn%128)*0.01,tofMatch1,tofMatch2,isTrue,nSigmaComb,nSigmaComb2,deltaphi1,deltaphi2,fPsi};
871 Float_t xTOfill2[] = {fPtLambdaC,KpTrack->Eta(),fPtKn,fPtKp,(fPidKn%128)*0.01,(fPidKp%128)*0.01,tofMatch2,tofMatch1,isTrue,nSigmaComb2,nSigmaComb,deltaphi2,deltaphi1,fPsi};
872
873
874 Int_t ipt = 0;
875 while(fPtLambdaMin[ipt] < fPtLambdaC && ipt < nPtBin){
876 ipt++;
877 }
878 ipt--;
879 if(ipt < 0) ipt = 0; // just to be sure
880
881 if(TMath::Abs(fEtaLambda) < 0.8 && fPtKp > 0.3 && fPtKn > 0.3){
882 if(isLambda) fContPid->Fill(0,fMassV0,fCentrality,xTOfill);
883 else fContPid2->Fill(0,fMassV0,fCentrality,xTOfill2);
520899fb 884
885
886
887 if(fPIDuserCut){
888 Float_t xUser[] = {KpTrack->Eta(),fPtKp,isTrue,0,deltaphi1,fPsi};
889
890 if(fPIDuserCut->IsSelected(KpTrack,AliPID::kPion)){ // to be filled
891 xUser[3] = 1;
892 } else if(fPIDuserCut->IsSelected(KpTrack,AliPID::kKaon)){
893 xUser[3] = 2;
894 } else if(fPIDuserCut->IsSelected(KpTrack,AliPID::kProton)){
895 xUser[3] = 3;
896 }
897
898 if(isLambda) fContUser->Fill(0,fMassV0,fCentrality,xUser);
899 else fContUser2->Fill(0,fMassV0,fCentrality,xUser);
900 }
a8ad4709 901 }
902
903
904 }
905 } // end analysi Lambda
906}
907
908//_____________________________________________________________________________
909Float_t AliAnalysisTaskLambdaBayes::GetVertex(AliAODEvent* aod) const
910{
911
912 Float_t zvtx = -999;
913
914 const AliAODVertex* vtxAOD = aod->GetPrimaryVertex();
915 if (!vtxAOD)
916 return zvtx;
917 if(vtxAOD->GetNContributors()>0)
918 zvtx = vtxAOD->GetZ();
919
920 return zvtx;
921}
922//_____________________________________________________________________________
923void AliAnalysisTaskLambdaBayes::Terminate(Option_t *)
924{
925 // Terminate loop
926 Printf("Terminate()");
927}
928//=======================================================================
929void AliAnalysisTaskLambdaBayes::SelectLambda(){
930 fNLambda=0;
931 fNpPos=0;
932 fNpNeg=0;
933
934 Int_t nV0s = fOutputAOD->GetNumberOfV0s();
935 AliAODv0 *myV0;
936 Double_t dMASS=0.0;
937 for (Int_t i=0; i!=nV0s; ++i) {
938 myV0 = (AliAODv0*) fOutputAOD->GetV0(i);
939 if(!myV0) continue;
940 if(myV0->Pt()<0.1 || TMath::Abs(myV0->Eta()) > 0.8) continue; // skipping low momentum
941 Int_t pass = PassesAODCuts(myV0,fOutputAOD,1); // check for Lambda
942 if(pass) {
943 if(pass==1) dMASS = myV0->MassLambda();
944 else dMASS = myV0->MassAntiLambda();
945
f4c4ac9e 946 Float_t massKs = myV0->MassK0Short();
947
948 if(TMath::Abs(dMASS-1.115)/0.005 < 8 && TMath::Abs(massKs - 0.497)/0.005 > 8){
a8ad4709 949 AliAODTrack *iT=(AliAODTrack*) myV0->GetDaughter(0); // positive
950 AliAODTrack *jT=(AliAODTrack*) myV0->GetDaughter(1); // negative
951 if(iT->Charge()<0){
952 iT=(AliAODTrack*) myV0->GetDaughter(1); // positive
953 jT=(AliAODTrack*) myV0->GetDaughter(0); // negative
954 }
955 fPhiLambda[fNLambda] = myV0->Phi();
956 fPtLambda[fNLambda] = myV0->Pt();
957 fIpP[fNLambda] = FindDaugheterIndex(iT);
958 fIpN[fNLambda] = FindDaugheterIndex(jT);
959
960 if(pass==2){
961 Int_t itemp = fIpP[fNLambda];
962 fIpP[fNLambda] = fIpN[fNLambda];
963 fIpN[fNLambda] = itemp;
964 }
965
966 fMassLambda[fNLambda] = dMASS;
967 if(fIpP[fNLambda] > -1 && fIpN[fNLambda] > -1){
968 fNLambda++;
969 }
970
971 }
972 }
973 }
974}
975
976//=======================================================================
977Int_t AliAnalysisTaskLambdaBayes::PassesAODCuts(AliAODv0 *myV0, AliAODEvent *tAOD,Int_t specie)
978{
979 if (myV0->GetOnFlyStatus() ) return 0;
980
981 //the following is needed in order to evualuate track-quality
982 AliAODTrack *iT, *jT;
983 AliAODVertex *vV0s = myV0->GetSecondaryVtx();
984 Double_t pos[3],cov[6];
985 vV0s->GetXYZ(pos);
986 vV0s->GetCovarianceMatrix(cov);
987 const AliESDVertex vESD(pos,cov,100.,100);
988
989 // TESTING CHARGE
990 int iPos, iNeg;
991 iT=(AliAODTrack*) myV0->GetDaughter(0);
992 if(iT->Charge()>0) {
993 iPos = 0; iNeg = 1;
994 } else {
995 iPos = 1; iNeg = 0;
996 }
997 // END OF TEST
998
999 iT=(AliAODTrack*) myV0->GetDaughter(iPos); // positive
1000
1001 jT=(AliAODTrack*) myV0->GetDaughter(iNeg); // negative
1002
1003 Bool_t trkFlag = iT->TestFilterBit(fFilterBit);
1004 if(!trkFlag) return 0;
1005 Bool_t trkFlag2 = jT->TestFilterBit(fFilterBit);
1006 if(!trkFlag2) return 0;
1007
1008 Double_t pvertex[3];
1009 pvertex[0]=tAOD->GetPrimaryVertex()->GetX();
1010 pvertex[1]=tAOD->GetPrimaryVertex()->GetY();
1011 pvertex[2]=tAOD->GetPrimaryVertex()->GetZ();
1012
1013 Double_t dDL=myV0->DecayLengthV0( pvertex );
1014 if(dDL < 0.5 || dDL > 25) return 0;
1015
1016 Double_t dDCA=myV0->DcaV0Daughters();
1017 if(dDCA >0.5) return 0;
1018
1019 Double_t dCTP=myV0->CosPointingAngle( pvertex );
1020 if(dCTP < -1) return 0;
1021
1022// AliESDtrack ieT( iT );
1023// Double_t dD0P=ieT.GetD(pvertex[0],pvertex[1],tAOD->GetMagneticField());
1024// if(TMath::Abs(dD0P) < 0]) return 0;
1025
1026// AliESDtrack jeT( jT );
1027// Double_t dD0M=jeT.GetD(pvertex[0],pvertex[1],tAOD->GetMagneticField());
1028// if(TMath::Abs(dD0M) < 0) return 0;
1029
1030// Double_t dD0D0=dD0P*dD0M;
1031// if(dD0D0>0) return 0;
1032
1033// Double_t dETA=myV0->Eta();
1034// if(dETA <-0.8) return 0;
1035// if(dETA >0.8) return 0;
1036
1037// Double_t dQT=myV0->PtArmV0();
1038// if(specie==0) if(dQT<???) return 0;
1039
1040 Double_t dALPHA=myV0->AlphaV0(); // AlphaV0 -> AODRecoDecat::Alpha -> return 1.-2./(1.+QlProng(0)/QlProng(1));
1041 if(myV0->ChargeProng(iPos)<0) dALPHA = -dALPHA; // protects for a change in convention
1042
1043 if(specie==1 && dALPHA<0) return 2; // antilambda
1044 return 1; // Ks or lambda
1045}
1046//-------------------------------------------------------------------------
1047Int_t AliAnalysisTaskLambdaBayes::FindDaugheterIndex(AliAODTrack *trk){
1048 Int_t ntrack = fOutputAOD->GetNumberOfTracks();
1049
1050 for(Int_t i=0;i < ntrack;i++){ // loop on tracks
1051 AliAODTrack* track = fOutputAOD->GetTrack(i);
1052 if(track == trk) return i;
1053 }
1054
1055 printf("Daughter for %p not found\n",trk);
1056 return -1;
1057}
1058//-------------------------------------------------------------------------
1059Int_t AliAnalysisTaskLambdaBayes::IsChannelValid(Float_t etaAbs){
1060 if(!fIsMC) return 1; // used only on MC
1061
9c668341 1062 if( fTypeCol==2){ // LHC10h or LHC11h because of TOF matching window at 3 cm
a8ad4709 1063 Int_t channel = Int_t(4334.09 - 4758.36 * etaAbs -1989.71 * etaAbs*etaAbs + 1957.62*etaAbs*etaAbs*etaAbs);
1064
1065 if(!(channel%20)) return 0; // 5% additional loss in MC
1066 }
1067
1068 return 1;
1069}