]> 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
It makes the addition in AOD of MC tracks associated to muons parametrizable from...
[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
788ca873 53 fProbTPCnSigma(),\r
54 fProbTOFnSigma(),\r
55 fProbTPCTOFnSigmaTPC(),\r
56 fProbTPC(),\r
57 fProbTOF(),\r
58 fProbTPCTOF(),\r
59 fPriors(),\r
60 fProbTPCTOFnSigTPCMom(),\r
61 fProbTPCnSigTPCMom(),\r
62 fProbTOFnSigTOFMom(),\r
63 fPriorsUsed(),\r
f5d6bdb0 64 fPIDResponse(0x0),\r
65 fPIDCombined(0x0),\r
66 fTrackCuts(0x0),\r
67 fTrackFilter(0x0)\r
68{\r
69 //\r
70 // Constructor\r
71 //\r
72}\r
73\r
74//_________________________________________________________________________________\r
75AliAnalysisTaskPIDCombined::AliAnalysisTaskPIDCombined(const char *name) :\r
76 AliAnalysisTaskSE(name),\r
77 fHistList(),\r
788ca873 78 fProbTPCnSigma(),\r
79 fProbTOFnSigma(),\r
80 fProbTPCTOFnSigmaTPC(),\r
81 fProbTPC(),\r
82 fProbTOF(),\r
83 fProbTPCTOF(),\r
84 fPriors(),\r
85 fProbTPCTOFnSigTPCMom(),\r
86 fProbTPCnSigTPCMom(),\r
87 fProbTOFnSigTOFMom(),\r
88 fPriorsUsed(),\r
f5d6bdb0 89 fPIDResponse(0x0),\r
90 fPIDCombined(0x0),\r
91 fTrackCuts(0x0),\r
92 fTrackFilter(0x0)\r
93{\r
94 //\r
95 // Constructor\r
96 //\r
97 DefineInput(0,TChain::Class());\r
98 DefineOutput(1, TList::Class());\r
99}\r
100\r
101//_________________________________________________________________________________\r
102void AliAnalysisTaskPIDCombined::UserCreateOutputObjects()\r
103{\r
104 //\r
105 // Initialise the framework objects\r
106 //\r
107\r
108\r
109 // ------- track cuts\r
110 fTrackCuts = new AliESDtrackCuts("fTrackCuts", "Standard");\r
111 fTrackCuts->SetAcceptKinkDaughters(kFALSE);\r
112 fTrackCuts->SetMinNClustersTPC(80);\r
113 fTrackCuts->SetMaxChi2PerClusterTPC(4);\r
114 fTrackCuts->SetMaxDCAToVertexXY(3);\r
115 fTrackCuts->SetMaxDCAToVertexZ(3);\r
116 fTrackCuts->SetRequireTPCRefit(kTRUE);\r
117 fTrackCuts->SetRequireITSRefit(kTRUE);\r
118 fTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);\r
119 fTrackFilter = new AliAnalysisFilter("trackFilter");\r
120 fTrackFilter->AddCuts(fTrackCuts);\r
121\r
122\r
123\r
124 // ------- setup PIDCombined\r
125 fPIDCombined=new AliPIDCombined;\r
126 fPIDCombined->SetDefaultTPCPriors();\r
127 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC+AliPIDResponse::kDetTOF);\r
128\r
129 // no light nuclei\r
130 fPIDCombined->SetSelectedSpecies(AliPID::kSPECIES);\r
131\r
132\r
133 fHistList.Add(new TH1D("nEvents","Number of Evnets;Selection",2,0,2));\r
134\r
135 for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec){\r
136\r
137 \r
138 fProbTPC[ispec]=new TH2D(Form("prob%s_mom_TPC",AliPID::ParticleName(ispec)),\r
139 Form("%s probability vs. momentum;momentum;probability",AliPID::ParticleName(ispec)),\r
140 100,0.,20.,50,0.,1.);\r
141 fHistList.Add(fProbTPC[ispec]);\r
142 fProbTPCnSigma[ispec]=new TH2D(Form("prob%s_nSigma_TPC",AliPID::ParticleName(ispec)),\r
143 Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
144 20,-5.,5.,50,0.,1.);\r
145 fHistList.Add(fProbTPCnSigma[ispec]);\r
146\r
147 for (Int_t ibin=0;ibin<kPtBins;ibin++) {\r
148 fProbTPCnSigTPCMom[ibin][ispec]=new TH2D(Form("prob%s_nSigma_TPC (%s)",AliPID::ParticleName(ispec),fgkBinMomDesc[ibin]),\r
149 Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
150 20,-5.,5.,50,0.,1.);\r
151 fHistList.Add(fProbTPCnSigTPCMom[ibin][ispec]);\r
152 }\r
153\r
154\r
155\r
156 fProbTOF[ispec]=new TH2D(Form("prob%s_mom_TOF",AliPID::ParticleName(ispec)),\r
157 Form("%s probability vs. momentum;momentum;probability",AliPID::ParticleName(ispec)),\r
158 100,0.,20.,50,0.,1.);\r
159 fHistList.Add(fProbTOF[ispec]);\r
160 fProbTOFnSigma[ispec]=new TH2D(Form("prob%s_nSigma_TOF",AliPID::ParticleName(ispec)),\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(fProbTOFnSigma[ispec]);\r
164 for (Int_t ibin=0;ibin<kPtBins;ibin++) {\r
165 fProbTOFnSigTOFMom[ibin][ispec]=new TH2D(Form("prob%s_nSigma_TOF (%s)",AliPID::ParticleName(ispec),fgkBinMomDesc[ibin]),\r
166 Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
167 20,-5.,5.,50,0.,1.);\r
168 fHistList.Add(fProbTOFnSigTOFMom[ibin][ispec]);\r
169 }\r
170\r
171\r
172\r
173 fProbTPCTOF[ispec]=new TH2D(Form("prob%s_mom_TPCTOF",AliPID::ParticleName(ispec)),\r
174 Form("%s probability vs. momentum;momentum;probability",AliPID::ParticleName(ispec)),\r
175 100,0.,20.,50,0.,1.);\r
176 fHistList.Add(fProbTPCTOF[ispec]);\r
177 fProbTPCTOFnSigmaTPC[ispec]=new TH2D(Form("probTPCTOF%s_nSigma_TPC",AliPID::ParticleName(ispec)),\r
178 Form("%s TPCTOF probability vs. n#sigmaTPC;n#sigma;probability",AliPID::ParticleName(ispec)),\r
179 20,-5.,5.,50,0.,1.);\r
180 fHistList.Add(fProbTPCTOFnSigmaTPC[ispec]);\r
181 for (Int_t ibin=0;ibin<kPtBins;ibin++) {\r
182 fProbTPCTOFnSigTPCMom[ibin][ispec]=new TH2D(Form("probTPCTOF%s_nSigma_TPC (%s)",AliPID::ParticleName(ispec),fgkBinMomDesc[ibin]),\r
183 Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
184 20,-5.,5.,50,0.,1.);\r
185 fHistList.Add(fProbTPCTOFnSigTPCMom[ibin][ispec]);\r
186 }\r
187\r
188\r
189\r
190 // basic priors\r
191 fPriors[ispec]=new TH1F(Form("%s_priors",AliPID::ParticleName(ispec)),\r
192 Form("%s priors vs momentum",AliPID::ParticleName(ispec)),\r
193 100,0.,20.);\r
194 fHistList.Add(fPriors[ispec]);\r
195 switch (ispec) {\r
196 case AliPID::kElectron:\r
197 for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.02);\r
198 break;\r
199 case AliPID::kMuon:\r
200 for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.02);\r
201 break;\r
202 case AliPID::kPion:\r
203 for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.56);\r
204 break;\r
205 case AliPID::kKaon:\r
206 for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.20);\r
207 break;\r
208 case AliPID::kProton:\r
209 for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.20);\r
210 break;\r
211 default:\r
212 break;\r
213 }\r
214 fPIDCombined->SetPriorDistribution((AliPID::EParticleType)ispec,fPriors[ispec]);\r
215\r
216 // priors used\r
217 fPriorsUsed[ispec] = new TH2D(Form("%s_priorsUsed",AliPID::ParticleName(ispec)),\r
218 Form("%s priors vs transverse momentum;p_{t} (GeV/c);priors",AliPID::ParticleName(ispec)),\r
219 100,0.,20.,101,0,1.01); \r
220 fHistList.Add(fPriorsUsed[ispec]);\r
221 }\r
222\r
223\r
224\r
225 fHistList.SetOwner();\r
226 PostData(1,&fHistList);\r
227\r
228\r
229}\r
230\r
231//_________________________________________________________________________________\r
232void AliAnalysisTaskPIDCombined::UserExec(Option_t *)\r
233{\r
234 //\r
235 // Main loop. Called for every event\r
236 //\r
237 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();\r
238 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());\r
239 fPIDResponse=inputHandler->GetPIDResponse();\r
240 if (!fPIDResponse) AliFatal("This Task needs the PID response attached to the inputHandler");\r
241\r
242 // Printf(" ---------------------- UserExec PID task ---------------------");\r
243 \r
244 FillHistogram("nEvents",0.);\r
245\r
246 AliVEvent *event=InputEvent();\r
247 AliVTrack *track=0x0;\r
248 Int_t ntracks=event->GetNumberOfTracks();\r
249\r
250 Double_t probTPC[AliPID::kSPECIES]={0.};\r
251 Double_t probTOF[AliPID::kSPECIES]={0.};\r
252 Double_t probTPCTOF[AliPID::kSPECIES]={0.};\r
253 \r
254 //loop over all tracks\r
255 for (Int_t itrack=0; itrack<ntracks; ++itrack){\r
256\r
257 track=(AliVTrack*)event->GetTrack(itrack);\r
258\r
259 if ( fTrackFilter->IsSelected(track) ) {\r
260\r
261 Double_t mom=track->GetTPCmomentum();\r
262 Int_t ibin=GetMomBin(mom);\r
263\r
264 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC);\r
265 UInt_t detUsed = fPIDCombined->ComputeProbabilities(track, fPIDResponse, probTPC);\r
266 \r
267 if (detUsed != 0) { // TPC is available\r
268 for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec) {\r
269 Double_t nSigmaTPC = fPIDResponse->NumberOfSigmasTPC(track,(AliPID::EParticleType)ispec);\r
270 fProbTPC[ispec]->Fill(mom,probTPC[ispec]);\r
271 fProbTPCnSigma[ispec]->Fill(nSigmaTPC,probTPC[ispec]);\r
272 fProbTPCnSigTPCMom[ibin][ispec]->Fill(nSigmaTPC,probTPC[ispec]);\r
273 }\r
274\r
275 // compute priors for TPC+TOF, even if we ask just TOF for PID\r
276 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTOF);\r
277 detUsed = fPIDCombined->ComputeProbabilities(track, fPIDResponse, probTOF);\r
278 Double_t priors[5]; // check priors used for TOF\r
279 fPIDCombined->GetPriors(track,priors,fPIDResponse,detUsed);\r
280 for(Int_t ispec=0;ispec<5;ispec++) fPriorsUsed[ispec]->Fill(TMath::Abs(track->Pt()),priors[ispec]);\r
281\r
282 if (detUsed != 0) { // TOF is available\r
283 for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec) {\r
284 Double_t nSigmaTOF = fPIDResponse->NumberOfSigmasTOF(track,(AliPID::EParticleType)ispec);\r
285 fProbTOF[ispec]->Fill(mom,probTOF[ispec]);\r
286 fProbTOFnSigma[ispec]->Fill(nSigmaTOF,probTOF[ispec]);\r
287 fProbTOFnSigTOFMom[ibin][ispec]->Fill(nSigmaTOF,probTOF[ispec]);\r
288 }\r
289 }\r
290\r
291 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTOF|AliPIDResponse::kDetTPC);\r
292 detUsed = fPIDCombined->ComputeProbabilities(track, fPIDResponse, probTPCTOF);\r
293 if (detUsed == (UInt_t)fPIDCombined->GetDetectorMask() ) {\r
294 for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec) {\r
295 Double_t nSigmaTPC = fPIDResponse->NumberOfSigmasTPC(track,(AliPID::EParticleType)ispec);\r
296 fProbTPCTOF[ispec]->Fill(mom,probTPCTOF[ispec]);\r
297 fProbTPCTOFnSigmaTPC[ispec]->Fill(nSigmaTPC,probTPCTOF[ispec]);\r
298 fProbTPCTOFnSigTPCMom[ibin][ispec]->Fill(nSigmaTPC,probTPCTOF[ispec]);\r
299 }\r
300 }\r
301\r
302 }\r
303\r
304\r
305 }\r
306 }\r
307\r
308 PostData(1, &fHistList);\r
309}\r
310\r
311//_________________________________________________________________________________\r
312void AliAnalysisTaskPIDCombined::FillHistogram(const char* name, Double_t x, Double_t weight)\r
313{\r
314 //\r
315 // Fill 1D histogram by name\r
316 //\r
317 ((TH1*)fHistList.FindObject(name))->Fill(x,weight);\r
318}\r
319\r
320//_________________________________________________________________________________\r
321void AliAnalysisTaskPIDCombined::FillHistogram(const char* name, Double_t x, Double_t y, Double_t weight)\r
322{\r
323 //\r
324 // Fill 2D histogram by name\r
325 //\r
326 ((TH2*)fHistList.FindObject(name))->Fill(x,y,weight);\r
327}\r
328\r
329\r
330//_________________________________________________________________________________\r
331Int_t AliAnalysisTaskPIDCombined::GetMomBin(Float_t mom)\r
332{\r
333 //\r
334 // Given momentum return histogram to be filled\r
335 //\r
336 if (mom>0. && mom < 0.5) return 0;\r
337 if (mom>=0.5 && mom < 0.7) return 1;\r
338 if (mom>=0.7 && mom < 1.0) return 2;\r
339 if (mom>=1.0 && mom < 1.5) return 3;\r
340 if (mom>=1.5 && mom < 2.0) return 4;\r
341 return kPtBins-1;\r
342}\r
343\r