]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ANALYSIS/AliAnalysisTaskPIDCombined.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Coverity 18637
[u/mrichter/AliRoot.git] / ANALYSIS / AliAnalysisTaskPIDCombined.cxx
CommitLineData
f5d6bdb0 1/*************************************************************************\r
2* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *\r
3* *\r
4* Author: The ALICE Off-line Project. *\r
5* Contributors are mentioned in the code where appropriate. *\r
6* *\r
7* Permission to use, copy, modify and distribute this software and its *\r
8* documentation strictly for non-commercial purposes is hereby granted *\r
9* without fee, provided that the above copyright notice appears in all *\r
10* copies and that both the copyright notice and this permission notice *\r
11* appear in the supporting documentation. The authors make no claims *\r
12* about the suitability of this software for any purpose. It is *\r
13* provided "as is" without express or implied warranty. *\r
14**************************************************************************/\r
15\r
16///////////////////////////////////////////////////////////////////////////\r
17// //\r
18// Basic Analysis Task //\r
19// for PID Analysis //\r
20// //\r
21///////////////////////////////////////////////////////////////////////////\r
22\r
23#include <TChain.h>\r
24#include <TH1D.h>\r
25#include <TH2D.h>\r
26\r
27#include <AliVEvent.h>\r
28#include <AliInputEventHandler.h>\r
29#include <AliAnalysisManager.h>\r
30\r
31#include <AliLog.h>\r
32#include <AliPID.h>\r
33#include <AliPIDCombined.h>\r
34#include <AliPIDResponse.h>\r
35\r
36#include "AliAnalysisTaskPIDCombined.h"\r
37\r
38const char *AliAnalysisTaskPIDCombined::fgkBinMomDesc[AliAnalysisTaskPIDCombined::kPtBins] = {\r
39 " 0 <= p < 0.5 GeV/c",\r
40 " 0.5 <= p < 0.7 GeV/c",\r
41 " 0.7 <= p < 1.0 GeV/c",\r
42 " 1.0 <= p < 1.5 GeV/c",\r
43 " 1.5 <= p < 2.0 GeV/c",\r
44 " p >= 2.0 GeV/c"\r
45};\r
46\r
47ClassImp(AliAnalysisTaskPIDCombined)\r
48\r
49//_________________________________________________________________________________\r
50AliAnalysisTaskPIDCombined::AliAnalysisTaskPIDCombined() :\r
51 AliAnalysisTaskSE(),\r
52 fHistList(),\r
53 fPIDResponse(0x0),\r
54 fPIDCombined(0x0),\r
55 fTrackCuts(0x0),\r
56 fTrackFilter(0x0)\r
57{\r
58 //\r
59 // Constructor\r
60 //\r
61}\r
62\r
63//_________________________________________________________________________________\r
64AliAnalysisTaskPIDCombined::AliAnalysisTaskPIDCombined(const char *name) :\r
65 AliAnalysisTaskSE(name),\r
66 fHistList(),\r
67 fPIDResponse(0x0),\r
68 fPIDCombined(0x0),\r
69 fTrackCuts(0x0),\r
70 fTrackFilter(0x0)\r
71{\r
72 //\r
73 // Constructor\r
74 //\r
75 DefineInput(0,TChain::Class());\r
76 DefineOutput(1, TList::Class());\r
77}\r
78\r
79//_________________________________________________________________________________\r
80void AliAnalysisTaskPIDCombined::UserCreateOutputObjects()\r
81{\r
82 //\r
83 // Initialise the framework objects\r
84 //\r
85\r
86\r
87 // ------- track cuts\r
88 fTrackCuts = new AliESDtrackCuts("fTrackCuts", "Standard");\r
89 fTrackCuts->SetAcceptKinkDaughters(kFALSE);\r
90 fTrackCuts->SetMinNClustersTPC(80);\r
91 fTrackCuts->SetMaxChi2PerClusterTPC(4);\r
92 fTrackCuts->SetMaxDCAToVertexXY(3);\r
93 fTrackCuts->SetMaxDCAToVertexZ(3);\r
94 fTrackCuts->SetRequireTPCRefit(kTRUE);\r
95 fTrackCuts->SetRequireITSRefit(kTRUE);\r
96 fTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);\r
97 fTrackFilter = new AliAnalysisFilter("trackFilter");\r
98 fTrackFilter->AddCuts(fTrackCuts);\r
99\r
100\r
101\r
102 // ------- setup PIDCombined\r
103 fPIDCombined=new AliPIDCombined;\r
104 fPIDCombined->SetDefaultTPCPriors();\r
105 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC+AliPIDResponse::kDetTOF);\r
106\r
107 // no light nuclei\r
108 fPIDCombined->SetSelectedSpecies(AliPID::kSPECIES);\r
109\r
110\r
111 fHistList.Add(new TH1D("nEvents","Number of Evnets;Selection",2,0,2));\r
112\r
113 for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec){\r
114\r
115 \r
116 fProbTPC[ispec]=new TH2D(Form("prob%s_mom_TPC",AliPID::ParticleName(ispec)),\r
117 Form("%s probability vs. momentum;momentum;probability",AliPID::ParticleName(ispec)),\r
118 100,0.,20.,50,0.,1.);\r
119 fHistList.Add(fProbTPC[ispec]);\r
120 fProbTPCnSigma[ispec]=new TH2D(Form("prob%s_nSigma_TPC",AliPID::ParticleName(ispec)),\r
121 Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
122 20,-5.,5.,50,0.,1.);\r
123 fHistList.Add(fProbTPCnSigma[ispec]);\r
124\r
125 for (Int_t ibin=0;ibin<kPtBins;ibin++) {\r
126 fProbTPCnSigTPCMom[ibin][ispec]=new TH2D(Form("prob%s_nSigma_TPC (%s)",AliPID::ParticleName(ispec),fgkBinMomDesc[ibin]),\r
127 Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
128 20,-5.,5.,50,0.,1.);\r
129 fHistList.Add(fProbTPCnSigTPCMom[ibin][ispec]);\r
130 }\r
131\r
132\r
133\r
134 fProbTOF[ispec]=new TH2D(Form("prob%s_mom_TOF",AliPID::ParticleName(ispec)),\r
135 Form("%s probability vs. momentum;momentum;probability",AliPID::ParticleName(ispec)),\r
136 100,0.,20.,50,0.,1.);\r
137 fHistList.Add(fProbTOF[ispec]);\r
138 fProbTOFnSigma[ispec]=new TH2D(Form("prob%s_nSigma_TOF",AliPID::ParticleName(ispec)),\r
139 Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
140 20,-5.,5.,50,0.,1.);\r
141 fHistList.Add(fProbTOFnSigma[ispec]);\r
142 for (Int_t ibin=0;ibin<kPtBins;ibin++) {\r
143 fProbTOFnSigTOFMom[ibin][ispec]=new TH2D(Form("prob%s_nSigma_TOF (%s)",AliPID::ParticleName(ispec),fgkBinMomDesc[ibin]),\r
144 Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
145 20,-5.,5.,50,0.,1.);\r
146 fHistList.Add(fProbTOFnSigTOFMom[ibin][ispec]);\r
147 }\r
148\r
149\r
150\r
151 fProbTPCTOF[ispec]=new TH2D(Form("prob%s_mom_TPCTOF",AliPID::ParticleName(ispec)),\r
152 Form("%s probability vs. momentum;momentum;probability",AliPID::ParticleName(ispec)),\r
153 100,0.,20.,50,0.,1.);\r
154 fHistList.Add(fProbTPCTOF[ispec]);\r
155 fProbTPCTOFnSigmaTPC[ispec]=new TH2D(Form("probTPCTOF%s_nSigma_TPC",AliPID::ParticleName(ispec)),\r
156 Form("%s TPCTOF probability vs. n#sigmaTPC;n#sigma;probability",AliPID::ParticleName(ispec)),\r
157 20,-5.,5.,50,0.,1.);\r
158 fHistList.Add(fProbTPCTOFnSigmaTPC[ispec]);\r
159 for (Int_t ibin=0;ibin<kPtBins;ibin++) {\r
160 fProbTPCTOFnSigTPCMom[ibin][ispec]=new TH2D(Form("probTPCTOF%s_nSigma_TPC (%s)",AliPID::ParticleName(ispec),fgkBinMomDesc[ibin]),\r
161 Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
162 20,-5.,5.,50,0.,1.);\r
163 fHistList.Add(fProbTPCTOFnSigTPCMom[ibin][ispec]);\r
164 }\r
165\r
166\r
167\r
168 // basic priors\r
169 fPriors[ispec]=new TH1F(Form("%s_priors",AliPID::ParticleName(ispec)),\r
170 Form("%s priors vs momentum",AliPID::ParticleName(ispec)),\r
171 100,0.,20.);\r
172 fHistList.Add(fPriors[ispec]);\r
173 switch (ispec) {\r
174 case AliPID::kElectron:\r
175 for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.02);\r
176 break;\r
177 case AliPID::kMuon:\r
178 for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.02);\r
179 break;\r
180 case AliPID::kPion:\r
181 for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.56);\r
182 break;\r
183 case AliPID::kKaon:\r
184 for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.20);\r
185 break;\r
186 case AliPID::kProton:\r
187 for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.20);\r
188 break;\r
189 default:\r
190 break;\r
191 }\r
192 fPIDCombined->SetPriorDistribution((AliPID::EParticleType)ispec,fPriors[ispec]);\r
193\r
194 // priors used\r
195 fPriorsUsed[ispec] = new TH2D(Form("%s_priorsUsed",AliPID::ParticleName(ispec)),\r
196 Form("%s priors vs transverse momentum;p_{t} (GeV/c);priors",AliPID::ParticleName(ispec)),\r
197 100,0.,20.,101,0,1.01); \r
198 fHistList.Add(fPriorsUsed[ispec]);\r
199 }\r
200\r
201\r
202\r
203 fHistList.SetOwner();\r
204 PostData(1,&fHistList);\r
205\r
206\r
207}\r
208\r
209//_________________________________________________________________________________\r
210void AliAnalysisTaskPIDCombined::UserExec(Option_t *)\r
211{\r
212 //\r
213 // Main loop. Called for every event\r
214 //\r
215 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();\r
216 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());\r
217 fPIDResponse=inputHandler->GetPIDResponse();\r
218 if (!fPIDResponse) AliFatal("This Task needs the PID response attached to the inputHandler");\r
219\r
220 // Printf(" ---------------------- UserExec PID task ---------------------");\r
221 \r
222 FillHistogram("nEvents",0.);\r
223\r
224 AliVEvent *event=InputEvent();\r
225 AliVTrack *track=0x0;\r
226 Int_t ntracks=event->GetNumberOfTracks();\r
227\r
228 Double_t probTPC[AliPID::kSPECIES]={0.};\r
229 Double_t probTOF[AliPID::kSPECIES]={0.};\r
230 Double_t probTPCTOF[AliPID::kSPECIES]={0.};\r
231 \r
232 //loop over all tracks\r
233 for (Int_t itrack=0; itrack<ntracks; ++itrack){\r
234\r
235 track=(AliVTrack*)event->GetTrack(itrack);\r
236\r
237 if ( fTrackFilter->IsSelected(track) ) {\r
238\r
239 Double_t mom=track->GetTPCmomentum();\r
240 Int_t ibin=GetMomBin(mom);\r
241\r
242 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC);\r
243 UInt_t detUsed = fPIDCombined->ComputeProbabilities(track, fPIDResponse, probTPC);\r
244 \r
245 if (detUsed != 0) { // TPC is available\r
246 for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec) {\r
247 Double_t nSigmaTPC = fPIDResponse->NumberOfSigmasTPC(track,(AliPID::EParticleType)ispec);\r
248 fProbTPC[ispec]->Fill(mom,probTPC[ispec]);\r
249 fProbTPCnSigma[ispec]->Fill(nSigmaTPC,probTPC[ispec]);\r
250 fProbTPCnSigTPCMom[ibin][ispec]->Fill(nSigmaTPC,probTPC[ispec]);\r
251 }\r
252\r
253 // compute priors for TPC+TOF, even if we ask just TOF for PID\r
254 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTOF);\r
255 detUsed = fPIDCombined->ComputeProbabilities(track, fPIDResponse, probTOF);\r
256 Double_t priors[5]; // check priors used for TOF\r
257 fPIDCombined->GetPriors(track,priors,fPIDResponse,detUsed);\r
258 for(Int_t ispec=0;ispec<5;ispec++) fPriorsUsed[ispec]->Fill(TMath::Abs(track->Pt()),priors[ispec]);\r
259\r
260 if (detUsed != 0) { // TOF is available\r
261 for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec) {\r
262 Double_t nSigmaTOF = fPIDResponse->NumberOfSigmasTOF(track,(AliPID::EParticleType)ispec);\r
263 fProbTOF[ispec]->Fill(mom,probTOF[ispec]);\r
264 fProbTOFnSigma[ispec]->Fill(nSigmaTOF,probTOF[ispec]);\r
265 fProbTOFnSigTOFMom[ibin][ispec]->Fill(nSigmaTOF,probTOF[ispec]);\r
266 }\r
267 }\r
268\r
269 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTOF|AliPIDResponse::kDetTPC);\r
270 detUsed = fPIDCombined->ComputeProbabilities(track, fPIDResponse, probTPCTOF);\r
271 if (detUsed == (UInt_t)fPIDCombined->GetDetectorMask() ) {\r
272 for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec) {\r
273 Double_t nSigmaTPC = fPIDResponse->NumberOfSigmasTPC(track,(AliPID::EParticleType)ispec);\r
274 fProbTPCTOF[ispec]->Fill(mom,probTPCTOF[ispec]);\r
275 fProbTPCTOFnSigmaTPC[ispec]->Fill(nSigmaTPC,probTPCTOF[ispec]);\r
276 fProbTPCTOFnSigTPCMom[ibin][ispec]->Fill(nSigmaTPC,probTPCTOF[ispec]);\r
277 }\r
278 }\r
279\r
280 }\r
281\r
282\r
283 }\r
284 }\r
285\r
286 PostData(1, &fHistList);\r
287}\r
288\r
289//_________________________________________________________________________________\r
290void AliAnalysisTaskPIDCombined::FillHistogram(const char* name, Double_t x, Double_t weight)\r
291{\r
292 //\r
293 // Fill 1D histogram by name\r
294 //\r
295 ((TH1*)fHistList.FindObject(name))->Fill(x,weight);\r
296}\r
297\r
298//_________________________________________________________________________________\r
299void AliAnalysisTaskPIDCombined::FillHistogram(const char* name, Double_t x, Double_t y, Double_t weight)\r
300{\r
301 //\r
302 // Fill 2D histogram by name\r
303 //\r
304 ((TH2*)fHistList.FindObject(name))->Fill(x,y,weight);\r
305}\r
306\r
307\r
308//_________________________________________________________________________________\r
309Int_t AliAnalysisTaskPIDCombined::GetMomBin(Float_t mom)\r
310{\r
311 //\r
312 // Given momentum return histogram to be filled\r
313 //\r
314 if (mom>0. && mom < 0.5) return 0;\r
315 if (mom>=0.5 && mom < 0.7) return 1;\r
316 if (mom>=0.7 && mom < 1.0) return 2;\r
317 if (mom>=1.0 && mom < 1.5) return 3;\r
318 if (mom>=1.5 && mom < 2.0) return 4;\r
319 return kPtBins-1;\r
320}\r
321\r