]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG1/ITS/AliAnalysisTaskITSsaTracks.cxx
add set and getter for neutral energy fraction
[u/mrichter/AliRoot.git] / PWG1 / ITS / AliAnalysisTaskITSsaTracks.cxx
CommitLineData
045571a6 1#include "AliAnalysisTaskSE.h"
2#include "AliTrackPointArray.h"
3#include "AliAnalysisManager.h"
4#include "AliAnalysisDataContainer.h"
5#include "AliESDEvent.h"
6#include "AliESDfriend.h"
7#include "AliStack.h"
8#include "AliPID.h"
9#include "AliITSPIDResponse.h"
10#include "AliMCEventHandler.h"
11#include "AliMCEvent.h"
12#include <TParticle.h>
13#include <TSystem.h>
14#include <TTree.h>
2afdea69 15#include <TNtuple.h>
045571a6 16#include <TH1F.h>
17#include <TH2F.h>
18#include <TChain.h>
19#include "AliESDInputHandlerRP.h"
20#include "AliAnalysisTaskITSsaTracks.h"
21
22/**************************************************************************
23 * Copyright(c) 1998-2012, ALICE Experiment at CERN, All rights reserved. *
24 * *
25 * Author: The ALICE Off-line Project. *
26 * Contributors are mentioned in the code where appropriate. *
27 * *
28 * Permission to use, copy, modify and distribute this software and its *
29 * documentation strictly for non-commercial purposes is hereby granted *
30 * without fee, provided that the above copyright notice appears in all *
31 * copies and that both the copyright notice and this permission notice *
32 * appear in the supporting documentation. The authors make no claims *
33 * about the suitability of this software for any purpose. It is *
34 * provided "as is" without express or implied warranty. *
35 **************************************************************************/
36
37//*************************************************************************
38// Implementation of class AliAnalysisTaskITSsaTracks
39// AliAnalysisTaskSE to extract QA and performance histos for ITS standalone tracks
40//
41//
42// Authors: L. Milano, milano@to.infn.it
43// F. Prino, prino@to.infn.it
44//
45//*************************************************************************
46
47ClassImp(AliAnalysisTaskITSsaTracks)
48//______________________________________________________________________________
49AliAnalysisTaskITSsaTracks::AliAnalysisTaskITSsaTracks() : AliAnalysisTaskSE("ITSsa resolution"),
50 fOutput(0),
51 fHistNEvents(0),
2afdea69 52 fHistMCPhiResid(0),
53 fHistPhiResid(0),
54 fNtupleTracks(0),
045571a6 55 fNPtBins(kMaxPtBins),
56 fMinITSpts(4),
57 fMinSPDpts(1),
2afdea69 58 fMinPtsforPid(3),
045571a6 59 fMinTPCpts(50),
2afdea69 60 fMaxITSChi2Clu(100.),
61 fFillNtuple(kFALSE),
045571a6 62 fReadMC(kFALSE),
63 fUseMCId(kFALSE)
64{
65 //
66 for(Int_t ilay=0; ilay<6;ilay++) fRequirePoint[ilay]=kFALSE;
67 DefineInput(0, TChain::Class());
68 DefineOutput(1, TList::Class());
69 const Int_t nbins = 29;
70 Double_t xbins[nbins+1]={0.06,0.08,0.10,0.12,0.14,0.16,0.18,0.20,0.25,0.30,
71 0.35,0.40,0.45,0.50,0.55,0.60,0.65,0.70,0.75,0.80,
72 0.85,0.90,0.95,1.00,1.20,1.40,1.60,1.80,1.90,2.00};
73 SetPtBins(nbins,xbins);
74}
75
76
77//___________________________________________________________________________
78AliAnalysisTaskITSsaTracks::~AliAnalysisTaskITSsaTracks(){
79 //
2afdea69 80 delete fHistNEvents;
81 for(Int_t iType=0; iType<kNtrackTypes; iType++){
82 delete fHistPt[iType];
83 delete fHistPtGood[iType];
84 delete fHistPtFake[iType];
85 delete fHistEtaPhi[iType];
86 delete fHistEtaPhiGood[iType];
87 delete fHistEtaPhiFake[iType];
519eed9b 88 delete fHistEtaPhiAny[iType];
89 delete fHistEtaPhi1SPD[iType];
90 delete fHistEtaPhi4Clu[iType];
91 delete fHistEtaPhi6Clu[iType];
2afdea69 92 delete fHistChi2[iType];
93 delete fHistChi2Good[iType];
94 delete fHistChi2Fake[iType];
95 delete fHistNclu[iType];
96 delete fHistNcluGood[iType];
97 delete fHistNcluFake[iType];
519eed9b 98 delete fHistdedxvsP2cls[iType];
99 delete fHistdedxvsP3cls[iType];
100 delete fHistdedxvsP4cls[iType];
2afdea69 101 }
102 for(Int_t iSpec=0; iSpec<kNspecies; iSpec++){
103 delete fHistPtTPCITS[iSpec];
104 delete fHistPtITSsa[iSpec];
105 delete fHistPtITSpureSA[iSpec];
106 delete fHistEtaPhiTPCITS[iSpec];
107 delete fHistEtaPhiITSsa[iSpec];
108 delete fHistEtaPhiITSpureSA[iSpec];
109 delete fHistNcluTPCITS[iSpec];
110 delete fHistNcluITSsa[iSpec];
111 delete fHistNcluITSpureSA[iSpec];
112 delete fHistd0rphiITSpureSA[iSpec];
113 delete fHistd0zITSpureSA[iSpec];
114 delete fHistCluInLayTPCITS[iSpec];
115 delete fHistCluInLayITSsa[iSpec];
116 delete fHistCluInLayITSpureSA[iSpec];
117 delete fHistOuterLayITSpureSA[iSpec];
118 delete fHistPtResid[iSpec];
119 delete fHistPtRelResid[iSpec];
120 delete fHistInvPtResid[iSpec];
121 delete fHistInvPtRelResid[iSpec];
122 delete fHistMCPtResid[iSpec];
123 delete fHistMCPtRelResid[iSpec];
124 delete fHistMCInvPtResid[iSpec];
125 delete fHistMCInvPtRelResid[iSpec];
126 }
127
128 delete fHistMCPhiResid;
129 delete fHistPhiResid;
130 delete fNtupleTracks;
045571a6 131 if (fOutput) {
132 delete fOutput;
133 fOutput = 0;
134 }
135}
136
137//________________________________________________________________________
138void AliAnalysisTaskITSsaTracks::SetPtBins(Int_t n, Double_t* lim){
139 // define pt bins for analysis
140 if(n>kMaxPtBins){
141 printf("Max. number of Pt bins = %d\n",kMaxPtBins);
142 return;
143 }else{
144 fNPtBins=n;
145 for(Int_t i=0; i<fNPtBins+1; i++) fPtLimits[i]=lim[i];
146 for(Int_t i=fNPtBins+1; i<kMaxPtBins+1; i++) fPtLimits[i]=99999999.;
147 }
148}
149//___________________________________________________________________________
150void AliAnalysisTaskITSsaTracks::UserCreateOutputObjects() {
151 // create output histos
152
153 fOutput = new TList();
154 fOutput->SetOwner();
155 fOutput->SetName("OutputHistos");
156
2afdea69 157 fNtupleTracks= new TNtuple("ntupleTracks","ntupleTracks","pt:px:py:pz:eta:phi:d0xy:d0z:dedx3:dedx4:dedx5:dedx6:truncmean:chi2:status:ITSrefit:TPCin:TPCrefit:ITSpureSA:nITSclu:nTPCclu:clumap:spdin:spdout:sddin:sddout:ssdin:ssdout:label:ptgen:pxgen:pygen:pzgen:etagen:phigen:ntracks");
158 // kinematics: pt, p eta,phi -> 6 variables
159 // impact parameters: d0xy, d0z -> 2 variables
160 // dE/dx: 4 Layers + trunc mean -> 5 variables
161 // chi2 and track status: -> 6 variables
162 // cluster infos: -> 9 variables
163 // MC info: -> 7 variables
164 // Multiplicity -> 1 variable
165 // Total: 36
166 fOutput->Add(fNtupleTracks);
167
045571a6 168 fHistNEvents = new TH1F("hNEvents", "Number of processed events",3,-0.5,2.5);
169 fHistNEvents->Sumw2();
170 fHistNEvents->SetMinimum(0);
171 fOutput->Add(fHistNEvents);
172
173 //binning for the dedx histogram
174 const Int_t hnbinsdedx=400;
175 Double_t hxmindedx = 0.01;
176 Double_t hxmaxdedx = 10;
177 Double_t hlogxmindedx = TMath::Log10(hxmindedx);
178 Double_t hlogxmaxdedx = TMath::Log10(hxmaxdedx);
179 Double_t hbinwidthdedx = (hlogxmaxdedx-hlogxmindedx)/hnbinsdedx;
180 Double_t hxbinsdedx[hnbinsdedx+1];
181 hxbinsdedx[0] = 0.01;
182 for (Int_t i=1;i<=hnbinsdedx;i++) {
183 hxbinsdedx[i] = hxmindedx + TMath::Power(10,hlogxmindedx+i*hbinwidthdedx);
184 }
185
2afdea69 186 TString tit[kNtrackTypes]={"TPCITS","ITSsa","ITSpureSA"};
187 for(Int_t iType=0; iType<kNtrackTypes; iType++){
188 fHistPt[iType]=new TH1F(Form("hPt%s",tit[iType].Data()),"",100,0.,2.);
189 fHistPt[iType]->Sumw2();
190 fOutput->Add(fHistPt[iType]);
191 fHistPtGood[iType]=new TH1F(Form("hPtGood%s",tit[iType].Data()),"",100,0.,2.);
192 fHistPtGood[iType]->Sumw2();
193 fOutput->Add(fHistPtGood[iType]);
194 fHistPtFake[iType]=new TH1F(Form("hPtFake%s",tit[iType].Data()),"",100,0.,2.);
195 fHistPtFake[iType]->Sumw2();
196 fOutput->Add(fHistPtFake[iType]);
197
198 fHistEtaPhi[iType] = new TH2F(Form("hEtaPhi%s",tit[iType].Data()),"",50,-1,1.,100,0.,2.*TMath::Pi());
199 fHistEtaPhi[iType]->Sumw2();
200 fOutput->Add(fHistEtaPhi[iType]);
201 fHistEtaPhiGood[iType] = new TH2F(Form("hEtaPhiGood%s",tit[iType].Data()),"",50,-1,1.,100,0.,2.*TMath::Pi());
202 fHistEtaPhiGood[iType]->Sumw2();
203 fOutput->Add(fHistEtaPhiGood[iType]);
204 fHistEtaPhiFake[iType] = new TH2F(Form("hEtaPhiFake%s",tit[iType].Data()),"",50,-1,1.,100,0.,2.*TMath::Pi());
205 fHistEtaPhiFake[iType]->Sumw2();
206 fOutput->Add(fHistEtaPhiFake[iType]);
207
519eed9b 208 fHistEtaPhiAny[iType] = new TH2F(Form("hEtaPhiAny%s",tit[iType].Data()),"",50,-1,1.,100,0.,2.*TMath::Pi());
209 fHistEtaPhiAny[iType]->Sumw2();
210 fOutput->Add(fHistEtaPhiAny[iType]);
211 fHistEtaPhi1SPD[iType] = new TH2F(Form("hEtaPhi1SPD%s",tit[iType].Data()),"",50,-1,1.,100,0.,2.*TMath::Pi());
212 fHistEtaPhi1SPD[iType]->Sumw2();
213 fOutput->Add(fHistEtaPhi1SPD[iType]);
214 fHistEtaPhi4Clu[iType] = new TH2F(Form("hEtaPhi4Clu%s",tit[iType].Data()),"",50,-1,1.,100,0.,2.*TMath::Pi());
215 fHistEtaPhi4Clu[iType]->Sumw2();
216 fOutput->Add(fHistEtaPhi4Clu[iType]);
217 fHistEtaPhi6Clu[iType] = new TH2F(Form("hEtaPhi6Clu%s",tit[iType].Data()),"",50,-1,1.,100,0.,2.*TMath::Pi());
218 fHistEtaPhi6Clu[iType]->Sumw2();
219 fOutput->Add(fHistEtaPhi6Clu[iType]);
220
2afdea69 221 fHistChi2[iType]=new TH1F(Form("hChi2%s",tit[iType].Data()),"",100,0.,10.);
222 fHistChi2[iType]->Sumw2();
223 fOutput->Add(fHistChi2[iType]);
224 fHistChi2Good[iType]=new TH1F(Form("hChi2Good%s",tit[iType].Data()),"",100,0.,10.);
225 fHistChi2Good[iType]->Sumw2();
226 fOutput->Add(fHistChi2Good[iType]);
227 fHistChi2Fake[iType]=new TH1F(Form("hChi2Fake%s",tit[iType].Data()),"",100,0.,10.);
228 fHistChi2Fake[iType]->Sumw2();
229 fOutput->Add(fHistChi2Fake[iType]);
230
231 fHistNclu[iType]=new TH1F(Form("hNclu%s",tit[iType].Data()),"",7,-0.5,6.5);
232 fHistNclu[iType]->Sumw2();
233 fOutput->Add(fHistNclu[iType]);
234 fHistNcluGood[iType]=new TH1F(Form("hNcluGood%s",tit[iType].Data()),"",7,-0.5,6.5);
235 fHistNcluGood[iType]->Sumw2();
236 fOutput->Add(fHistNcluGood[iType]);
237 fHistNcluFake[iType]=new TH1F(Form("hNcluFake%s",tit[iType].Data()),"",7,-0.5,6.5);
238 fHistNcluFake[iType]->Sumw2();
239 fOutput->Add(fHistNcluFake[iType]);
240
519eed9b 241 fHistdedxvsP2cls[iType] = new TH2F(Form("hdedxvsP2cls%s",tit[iType].Data()),"",hnbinsdedx,hxbinsdedx,900,0,1000);
242 fHistdedxvsP2cls[iType]->Sumw2();
243 fOutput->Add(fHistdedxvsP2cls[iType]);
244
245 fHistdedxvsP3cls[iType] = new TH2F(Form("hdedxvsP3cls%s",tit[iType].Data()),"",hnbinsdedx,hxbinsdedx,900,0,1000);
246 fHistdedxvsP3cls[iType]->Sumw2();
247 fOutput->Add(fHistdedxvsP3cls[iType]);
2afdea69 248
519eed9b 249 fHistdedxvsP4cls[iType] = new TH2F(Form("hdedxvsP4cls%s",tit[iType].Data()),"",hnbinsdedx,hxbinsdedx,900,0,1000);
250 fHistdedxvsP4cls[iType]->Sumw2();
251 fOutput->Add(fHistdedxvsP4cls[iType]);
2afdea69 252 }
253
254 //-----------------------------------------------------------
255
045571a6 256 TString spname[3]={"Pion","Kaon","Proton"};
257 TString hisname;
258 const Int_t nbins = fNPtBins;
259 Double_t xbins[nbins+1];
260 for(Int_t ibin=0; ibin<=nbins; ibin++) xbins[ibin]=fPtLimits[ibin];
261
045571a6 262 for(Int_t iSpec=0; iSpec<kNspecies; iSpec++){
263
264 hisname.Form("hPtTPCITS%s",spname[iSpec].Data());
265 fHistPtTPCITS[iSpec] = new TH1F(hisname.Data(),"",100,0.,2.);
266 fHistPtTPCITS[iSpec]->Sumw2();
267 fOutput->Add(fHistPtTPCITS[iSpec]);
268
269 hisname.Form("hPtITSsa%s",spname[iSpec].Data());
270 fHistPtITSsa[iSpec] = new TH1F(hisname.Data(),"",100,0.,2.);
271 fHistPtITSsa[iSpec]->Sumw2();
272 fOutput->Add(fHistPtITSsa[iSpec]);
273
274 hisname.Form("hPtITSpureSA%s",spname[iSpec].Data());
275 fHistPtITSpureSA[iSpec] = new TH1F(hisname.Data(),"",100,0.,2.);
276 fHistPtITSpureSA[iSpec]->Sumw2();
277 fOutput->Add(fHistPtITSpureSA[iSpec]);
278
279 //---
280
281 hisname.Form("hEtaPhiTPCITS%s",spname[iSpec].Data());
282 fHistEtaPhiTPCITS[iSpec] = new TH2F(hisname.Data(),"",50,-1,1.,50,0.,2.*TMath::Pi());
283 fHistEtaPhiTPCITS[iSpec]->Sumw2();
284 fOutput->Add(fHistEtaPhiTPCITS[iSpec]);
285
286 hisname.Form("hEtaPhiITSsa%s",spname[iSpec].Data());
287 fHistEtaPhiITSsa[iSpec] = new TH2F(hisname.Data(),"",50,-1,1.,50,0.,2.*TMath::Pi());
288 fHistEtaPhiITSsa[iSpec]->Sumw2();
289 fOutput->Add(fHistEtaPhiITSsa[iSpec]);
290
291 hisname.Form("hEtaPhiITSpureSA%s",spname[iSpec].Data());
292 fHistEtaPhiITSpureSA[iSpec] = new TH2F(hisname.Data(),"",50,-1,1.,50,0.,2.*TMath::Pi());
293 fHistEtaPhiITSpureSA[iSpec]->Sumw2();
294 fOutput->Add(fHistEtaPhiITSpureSA[iSpec]);
295
296 //---
297
298 hisname.Form("hNcluTPCITS%s",spname[iSpec].Data());
299 fHistNcluTPCITS[iSpec] = new TH2F(hisname.Data(),"",100,0.,2.,7,-0.5,6.5);
300 fHistNcluTPCITS[iSpec]->Sumw2();
301 fOutput->Add(fHistNcluTPCITS[iSpec]);
302
303 hisname.Form("hNcluITSsa%s",spname[iSpec].Data());
304 fHistNcluITSsa[iSpec] = new TH2F(hisname.Data(),"",100,0.,2.,7,-0.5,6.5);
305 fHistNcluITSsa[iSpec]->Sumw2();
306 fOutput->Add(fHistNcluITSsa[iSpec]);
307
308 hisname.Form("hNcluITSpureSA%s",spname[iSpec].Data());
309 fHistNcluITSpureSA[iSpec] = new TH2F(hisname.Data(),"",100,0.,2.,7,-0.5,6.5);
310 fHistNcluITSpureSA[iSpec]->Sumw2();
311 fOutput->Add(fHistNcluITSpureSA[iSpec]);
312
313 //---
314
315 hisname.Form("hd0rphiITSpureSA%s",spname[iSpec].Data());
316 fHistd0rphiITSpureSA[iSpec] = new TH2F(hisname.Data(),"",nbins,xbins,2000,-1,1);
317 fHistd0rphiITSpureSA[iSpec]->Sumw2();
318 fOutput->Add(fHistd0rphiITSpureSA[iSpec]);
319
320 hisname.Form("hd0zITSpureSA%s",spname[iSpec].Data());
321 fHistd0zITSpureSA[iSpec] = new TH2F(hisname.Data(),"",nbins,xbins,2000,-1,1);
322 fHistd0zITSpureSA[iSpec]->Sumw2();
323 fOutput->Add(fHistd0zITSpureSA[iSpec]);
324
325 //---
326
327 hisname.Form("hCluInLayTPCITS%s",spname[iSpec].Data());
328 fHistCluInLayTPCITS[iSpec] = new TH2F(hisname.Data(),"",100,0.,2.,7,-1.5,5.5);
329 fHistCluInLayTPCITS[iSpec]->Sumw2();
330 fOutput->Add(fHistCluInLayTPCITS[iSpec]);
331
332 hisname.Form("hCluInLayITSsa%s",spname[iSpec].Data());
333 fHistCluInLayITSsa[iSpec] = new TH2F(hisname.Data(),"",100,0.,2.,7,-1.5,5.5);
334 fHistCluInLayITSsa[iSpec]->Sumw2();
335 fOutput->Add(fHistCluInLayITSsa[iSpec]);
336
337 hisname.Form("hCluInLayITSpureSA%s",spname[iSpec].Data());
338 fHistCluInLayITSpureSA[iSpec] = new TH2F(hisname.Data(),"",100,0.,2.,7,-1.5,5.5);
339 fHistCluInLayITSpureSA[iSpec]->Sumw2();
340 fOutput->Add(fHistCluInLayITSpureSA[iSpec]);
341
342 hisname.Form("hOuterLayITSpureSA%s",spname[iSpec].Data());
343 fHistOuterLayITSpureSA[iSpec] = new TH2F(hisname.Data(),"",100,0.,2.,7,-1.5,5.5);
344 fHistOuterLayITSpureSA[iSpec]->Sumw2();
345 fOutput->Add(fHistOuterLayITSpureSA[iSpec]);
346
347 //---
348
349 hisname.Form("hPtResid%s",spname[iSpec].Data());
350 fHistPtResid[iSpec]=new TH2F(hisname.Data(),hisname.Data(),nbins,xbins,100,-1.,1.);
351 fHistPtResid[iSpec]->Sumw2();
352 fOutput->Add(fHistPtResid[iSpec]);
353 hisname.Form("hPtRelResid%s",spname[iSpec].Data());
354 fHistPtRelResid[iSpec]=new TH2F(hisname.Data(),hisname.Data(),nbins,xbins,100,-0.5,0.5);
355 fHistPtRelResid[iSpec]->Sumw2();
356 fOutput->Add(fHistPtRelResid[iSpec]);
357
358 hisname.Form("hInvPtResid%s",spname[iSpec].Data());
359 fHistInvPtResid[iSpec]=new TH2F(hisname.Data(),hisname.Data(),nbins,xbins,100,-1.,1.);
360 fHistInvPtResid[iSpec]->Sumw2();
361 fOutput->Add(fHistInvPtResid[iSpec]);
362 hisname.Form("hInvPtRelResid%s",spname[iSpec].Data());
363 fHistInvPtRelResid[iSpec]=new TH2F(hisname.Data(),hisname.Data(),nbins,xbins,100,-0.5,0.5);
364 fHistInvPtRelResid[iSpec]->Sumw2();
365 fOutput->Add(fHistInvPtRelResid[iSpec]);
366
367 hisname.Form("hMCPtResid%s",spname[iSpec].Data());
368 fHistMCPtResid[iSpec]=new TH2F(hisname.Data(),hisname.Data(),nbins,xbins,100,-1.,1.);
369 fHistMCPtResid[iSpec]->Sumw2();
370 fOutput->Add(fHistMCPtResid[iSpec]);
371 hisname.Form("hMCPtRelResid%s",spname[iSpec].Data());
372 fHistMCPtRelResid[iSpec]=new TH2F(hisname.Data(),hisname.Data(),nbins,xbins,100,-0.5,0.5);
373 fHistMCPtRelResid[iSpec]->Sumw2();
374 fOutput->Add(fHistMCPtRelResid[iSpec]);
375
376 hisname.Form("hMCInvPtResid%s",spname[iSpec].Data());
377 fHistMCInvPtResid[iSpec]=new TH2F(hisname.Data(),hisname.Data(),nbins,xbins,100,-1.,1.);
378 fHistMCInvPtResid[iSpec]->Sumw2();
379 fOutput->Add(fHistMCInvPtResid[iSpec]);
380 hisname.Form("hMCInvPtRelResid%s",spname[iSpec].Data());
381 fHistMCInvPtRelResid[iSpec]=new TH2F(hisname.Data(),hisname.Data(),nbins,xbins,100,-0.5,0.5);
382 fHistMCInvPtRelResid[iSpec]->Sumw2();
383 fOutput->Add(fHistMCInvPtRelResid[iSpec]);
384
385 }
386
2afdea69 387 fHistMCPhiResid=new TH2F("hMCPhiResid","",nbins,xbins,100,-0.5,0.5);
388 fHistMCPhiResid->Sumw2();
389 fOutput->Add(fHistMCPhiResid);
390 fHistPhiResid=new TH2F("hPhiResid","",nbins,xbins,100,-0.5,0.5);
391 fHistPhiResid->Sumw2();
392 fOutput->Add(fHistPhiResid);
393
394 PostData(1,fOutput);
395
045571a6 396}
397//______________________________________________________________________________
398void AliAnalysisTaskITSsaTracks::UserExec(Option_t *)
399{
400 //
401
402 AliESDEvent *esd = (AliESDEvent*) (InputEvent());
403
404
405 if(!esd) {
406 printf("AliAnalysisTaskSDDRP::Exec(): bad ESD\n");
407 return;
408 }
409
410
411 if(!ESDfriend()) {
412 printf("AliAnalysisTaskSDDRP::Exec(): bad ESDfriend\n");
413 return;
414 }
415
416 AliStack* stack=0;
417
418 if(fReadMC){
419 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
420 if (!eventHandler) {
421 Printf("ERROR: Could not retrieve MC event handler");
422 return;
423 }
424 AliMCEvent* mcEvent = eventHandler->MCEvent();
425 if (!mcEvent) {
426 Printf("ERROR: Could not retrieve MC event");
427 return;
428 }
429 stack = mcEvent->Stack();
430 if (!stack) {
431 Printf("ERROR: stack not available");
432 return;
433 }
434 }
435
045571a6 436
437 fHistNEvents->Fill(0);
438 const AliESDVertex *spdv=esd->GetPrimaryVertexSPD();
439 if(spdv->GetNContributors()<=0) return;
440 fHistNEvents->Fill(1);
2afdea69 441
442 const Int_t ntSize=36;
443 Float_t xnt[ntSize];
444
045571a6 445 Int_t ntracks = esd->GetNumberOfTracks();
446 for (Int_t iTrack=0; iTrack < ntracks; iTrack++) {
447 AliESDtrack * track = esd->GetTrack(iTrack);
448 if (!track) continue;
2afdea69 449 Float_t pttrack=track->Pt();
450 Float_t ptrack=track->P();
451 Float_t pxtrack=track->Px();
452 Float_t pytrack=track->Py();
453 Float_t pztrack=track->Pz();
454 Float_t impactXY=-999, impactZ=-999;
455 track->GetImpactParameters(impactXY, impactZ);
456 Double_t dedxlay[4];
457 track->GetITSdEdxSamples(dedxlay);
458 Float_t dedx=track->GetITSsignal();
459 Float_t chi2=track->GetITSchi2();
045571a6 460 Int_t status=track->GetStatus();
2afdea69 461 Bool_t isITSrefit=kFALSE;
462 Bool_t isTPCin=kFALSE;
463 Bool_t isTPCrefit=kFALSE;
464 Bool_t isPureSA=kFALSE;
465 if(status&AliESDtrack::kITSrefit) isITSrefit=kTRUE;
466 if(status&AliESDtrack::kTPCin) isTPCin=kTRUE;
467 if(status&AliESDtrack::kTPCrefit) isTPCrefit=kTRUE;
468 if(status&AliESDtrack::kITSpureSA) isPureSA=kTRUE;
045571a6 469 Bool_t isSA=kTRUE;
470 if(status&AliESDtrack::kTPCin) isSA=kFALSE;
471 Int_t nTPCclus=track->GetNcls(1);
472 Int_t nITSclus=track->GetNcls(0);
2afdea69 473 UChar_t clumap=track->GetITSClusterMap();
474 Int_t idet;
475 Float_t xloc,zloc;
476 Int_t statusLay[6];
477 for(Int_t iLay=0; iLay<6;iLay++) track->GetITSModuleIndexInfo(iLay,idet,statusLay[iLay],xloc,zloc);
478 Int_t trlabel=track->GetLabel();
479 Float_t ptgen=-999.;
480 Float_t pgen=-999.;
481 Float_t pxgen=-999.;
482 Float_t pygen=-999.;
483 Float_t pzgen=-999.;
484 Float_t etagen=-999.;
485 Float_t phigen=-999.;
486 if(fReadMC){
487 TParticle* part = stack->Particle(TMath::Abs(trlabel));
488 ptgen=part->Pt();
489 pgen=part->P();
490 pxgen=part->Px();
491 pygen=part->Py();
492 pzgen=part->Pz();
493 etagen=part->Eta();
494 phigen=part->Phi();
495 }
496
497 Int_t indexn=0;
498 xnt[indexn++]=pttrack;
499 xnt[indexn++]=pxtrack;
500 xnt[indexn++]=pytrack;
501 xnt[indexn++]=pztrack;
502 xnt[indexn++]=track->Eta();
503 xnt[indexn++]=track->Phi();
504 xnt[indexn++]=impactXY;
505 xnt[indexn++]=impactZ;
506 for(Int_t iLay=0; iLay<4; iLay++) xnt[indexn++]=dedxlay[iLay];
507 xnt[indexn++]=dedx;
508 xnt[indexn++]=chi2;
509 xnt[indexn++]=status;
510 xnt[indexn++]=isITSrefit;
511 xnt[indexn++]=isTPCin;
512 xnt[indexn++]=isTPCrefit;
513 xnt[indexn++]=isPureSA;
514 xnt[indexn++]=nITSclus;
515 xnt[indexn++]=nTPCclus;
516 xnt[indexn++]=clumap;
517 for(Int_t iLay=0; iLay<6; iLay++) xnt[indexn++]=(Float_t)statusLay[iLay];
518 xnt[indexn++]=trlabel;
519 xnt[indexn++]=ptgen;
520 xnt[indexn++]=pxgen;
521 xnt[indexn++]=pygen;
522 xnt[indexn++]=pzgen;
523 xnt[indexn++]=etagen;
524 xnt[indexn++]=phigen;
525 xnt[indexn++]=ntracks;
526
527 if(indexn>ntSize) printf("AliAnalysisTaskITSsaTracks: ERROR ntuple insexout of range\n");
528 if(fFillNtuple) fNtupleTracks->Fill(xnt);
529
530 if(!(status&AliESDtrack::kITSrefit)) continue;
519eed9b 531
532 Int_t iTrackType=-1;
533 if(status&AliESDtrack::kTPCin){
534 iTrackType=kTypeTPCITS;
535 }else{
536 if(status&AliESDtrack::kITSpureSA) iTrackType=kTypeITSpureSA;
537 else iTrackType=kTypeITSsa;
538 }
539 if(iTrackType<0 || iTrackType>=kNtrackTypes) continue;
045571a6 540
045571a6 541 Int_t nSPDPoints=0;
542 for(Int_t i=0; i<2; i++){
543 if(clumap&(1<<i)) ++nSPDPoints;
544 }
045571a6 545 Int_t nPointsForPid=0;
546 for(Int_t i=2; i<6; i++){
547 if(clumap&(1<<i)) ++nPointsForPid;
548 }
2afdea69 549
519eed9b 550 fHistEtaPhiAny[iTrackType]->Fill(track->Eta(),track->Phi());
551 if(nSPDPoints>=1) fHistEtaPhi1SPD[iTrackType]->Fill(track->Eta(),track->Phi());
552 if(nSPDPoints>=1 && nPointsForPid>=3) fHistEtaPhi4Clu[iTrackType]->Fill(track->Eta(),track->Phi());
553 if(nITSclus==6) fHistEtaPhi6Clu[iTrackType]->Fill(track->Eta(),track->Phi());
554
555
556 if(nITSclus<fMinITSpts)continue;
557 if((chi2/nITSclus) > fMaxITSChi2Clu) continue;
558 if(nSPDPoints<fMinSPDpts) continue;
559
2afdea69 560 fHistPt[iTrackType]->Fill(pttrack);
561 fHistEtaPhi[iTrackType]->Fill(track->Eta(),track->Phi());
519eed9b 562
2afdea69 563 fHistChi2[iTrackType]->Fill(chi2/nITSclus);
564 fHistNclu[iTrackType]->Fill(nITSclus);
519eed9b 565 if(nPointsForPid==2){
566 fHistdedxvsP2cls[iTrackType]->Fill(ptrack,dedx);
567 }
568 else if(nPointsForPid==3){
569 fHistdedxvsP3cls[iTrackType]->Fill(ptrack,dedx);
2afdea69 570 }
519eed9b 571 else if(nPointsForPid==4){
572 fHistdedxvsP4cls[iTrackType]->Fill(ptrack,dedx);
2afdea69 573 }
574 if(fReadMC){
575 if(trlabel<0){
576 fHistPtFake[iTrackType]->Fill(pttrack);
577 fHistEtaPhiFake[iTrackType]->Fill(track->Eta(),track->Phi());
578 fHistChi2Fake[iTrackType]->Fill(chi2/nITSclus);
579 fHistNcluFake[iTrackType]->Fill(nITSclus);
045571a6 580 }else{
2afdea69 581 fHistPtGood[iTrackType]->Fill(pttrack);
582 fHistEtaPhiGood[iTrackType]->Fill(track->Eta(),track->Phi());
583 fHistChi2Good[iTrackType]->Fill(chi2/nITSclus);
584 fHistNcluGood[iTrackType]->Fill(nITSclus);
045571a6 585 }
586 }
587
2afdea69 588 Int_t hadronSpecie=-1;
045571a6 589 if(fReadMC && fUseMCId){
045571a6 590 if(trlabel>=0){
2afdea69 591 TParticle* part = stack->Particle(trlabel);
592 Int_t pdg=TMath::Abs(part->GetPdgCode());
593 if(pdg==211) hadronSpecie=kPion;
594 else if(pdg==321) hadronSpecie=kKaon;
595 else if(pdg==2212) hadronSpecie=kProton;
045571a6 596 }
045571a6 597 }else{
2afdea69 598 if(nPointsForPid<fMinPtsforPid)continue;
045571a6 599 AliITSPIDResponse pidits(fReadMC);
600 Float_t nSigmaPion=pidits.GetNumberOfSigmas(ptrack,dedx,AliPID::kPion,nPointsForPid,isSA);
601 if(nSigmaPion>-2. && nSigmaPion<2.){
602 hadronSpecie=kPion;
603 }else{
604 Float_t nSigmaKaon=pidits.GetNumberOfSigmas(ptrack,dedx,AliPID::kKaon,nPointsForPid,isSA);
605 if(nSigmaKaon>-2. && nSigmaKaon<2.){
606 hadronSpecie=kKaon;
607 }else{
608 Float_t nSigmaProton=pidits.GetNumberOfSigmas(ptrack,dedx,AliPID::kProton,nPointsForPid,isSA);
609 if(nSigmaProton>-2. && nSigmaProton<2.){
610 hadronSpecie=kProton;
611 }
612 }
613 }
614 }
615 if(hadronSpecie<0) continue;
2afdea69 616 if(iTrackType==kTypeTPCITS){ // TPC+ITS tracks
045571a6 617 fHistPtTPCITS[hadronSpecie]->Fill(pttrack);
618 fHistEtaPhiTPCITS[hadronSpecie]->Fill(track->Eta(),track->Phi());
619 fHistNcluTPCITS[hadronSpecie]->Fill(pttrack,nITSclus);
620 fHistCluInLayTPCITS[hadronSpecie]->Fill(-1.);
621 for(Int_t iBit=0; iBit<6; iBit++){
622 if(clumap&(1<<iBit)) fHistCluInLayTPCITS[hadronSpecie]->Fill(pttrack,iBit);
623 }
2afdea69 624 }else if(iTrackType==kTypeITSpureSA){ // TPC+ITS tracks
625 fHistPtITSpureSA[hadronSpecie]->Fill(pttrack);
626 fHistEtaPhiITSpureSA[hadronSpecie]->Fill(track->Eta(),track->Phi());
627 fHistNcluITSpureSA[hadronSpecie]->Fill(pttrack,nITSclus);
628 fHistd0rphiITSpureSA[hadronSpecie]->Fill(pttrack,impactXY);
629 fHistd0zITSpureSA[hadronSpecie]->Fill(pttrack,impactZ);
630 fHistCluInLayITSpureSA[hadronSpecie]->Fill(-1.);
631 Int_t outerLay=-1;
632 for(Int_t iBit=0; iBit<6; iBit++){
633 if(clumap&(1<<iBit)){
634 fHistCluInLayITSpureSA[hadronSpecie]->Fill(pttrack,iBit);
635 if(iBit>outerLay) outerLay=iBit;
045571a6 636 }
045571a6 637 }
2afdea69 638 fHistOuterLayITSpureSA[hadronSpecie]->Fill(pttrack,outerLay);
639
045571a6 640
2afdea69 641 if(fReadMC){
642 if(trlabel>=0){
045571a6 643 TParticle* part = stack->Particle(trlabel);
2afdea69 644 ptgen=part->Pt();
045571a6 645 Float_t invpttrack=track->OneOverPt();
646 Float_t invptgen=0.;
647 if(ptgen>0.) invptgen=1./ptgen;
648 fHistMCPtResid[hadronSpecie]->Fill(pttrack,pttrack-ptgen);
649 fHistMCPtRelResid[hadronSpecie]->Fill(pttrack,(pttrack-ptgen)/ptgen);
650 fHistMCInvPtResid[hadronSpecie]->Fill(pttrack,invpttrack-invptgen);
651 fHistMCInvPtRelResid[hadronSpecie]->Fill(pttrack,(invpttrack-invptgen)/invptgen);
2afdea69 652 Float_t deltaphi=track->Phi()-part->Phi();
653 if(deltaphi<-TMath::Pi()) deltaphi+=2*TMath::Pi();
654 if(deltaphi>TMath::Pi()) deltaphi-=2*TMath::Pi();
655 fHistMCPhiResid->Fill(pttrack,deltaphi);
045571a6 656 }
2afdea69 657 }
658 }else if(iTrackType==kTypeITSsa){ // TPC+ITS tracks
659 fHistPtITSsa[hadronSpecie]->Fill(pttrack);
660 fHistEtaPhiITSsa[hadronSpecie]->Fill(track->Eta(),track->Phi());
661 fHistNcluITSsa[hadronSpecie]->Fill(pttrack,nITSclus);
662 fHistCluInLayITSsa[hadronSpecie]->Fill(-1.);
663 for(Int_t iBit=0; iBit<6; iBit++){
664 if(clumap&(1<<iBit)) fHistCluInLayITSsa[hadronSpecie]->Fill(pttrack,iBit);
045571a6 665 }
666 }
2afdea69 667
668 // match TPCITS with ITSpureSA
669 if(nITSclus<fMinITSpts || nTPCclus<fMinTPCpts) continue;
045571a6 670 Int_t idxMI[12],idxSA[12];
671 for(Int_t icl=0; icl<12; icl++){
672 idxMI[icl]=-1;
673 idxSA[icl]=-1;
674 }
675 Int_t ncls=track->GetClusters(0,idxMI);
676 if(fMinITSpts<6){
677 Bool_t accept=kTRUE;
678 for(Int_t ilay=0; ilay<6; ilay++){
679 if(fRequirePoint[ilay]){
680 Int_t mask = 1<<ilay;
681 if(!(clumap & mask)){
682 accept=kFALSE;
683 break;
684 }
685 }
686 }
687 if(!accept) continue;
688 }
689 // Sort
690 for(Int_t i=0;i<12;i++){
691 for(Int_t j=i+1;j<12;j++){
692 if(idxMI[j]>idxMI[i]){
693 Int_t tmp=idxMI[j];
694 idxMI[j]=idxMI[i];
695 idxMI[i]=tmp;
696 }
697 }
698 }
699 // for(Int_t i=0; i<12; i++) printf("%d ",idxMI[i]);
700 // printf("\n");
701 if(idxMI[0]<0 && idxMI[0]==idxMI[1]) continue;
702 Bool_t matched=kFALSE;
703 for (Int_t iTrack2 = 0; iTrack2 < ntracks; iTrack2++) {
704 if(matched) break;
705 if(iTrack2==iTrack) continue;
706 AliESDtrack* track2 = esd->GetTrack(iTrack2);
707 Int_t status2=track2->GetStatus();
708 if(!(status2&AliESDtrack::kITSrefit)) continue;
709 if(!(status2&AliESDtrack::kITSpureSA)) continue;
710 Int_t clumap2=track2->GetITSClusterMap();
711 Int_t nITSclus2=track2->GetNcls(0);
712 Int_t nTPCclus2=track2->GetNcls(1);
713 if(nITSclus2<fMinITSpts || nTPCclus2>0) continue;
714 Int_t ncls2=track2->GetClusters(0,idxSA);
715 if(ncls2!=ncls) continue;
716 if(fMinITSpts<6){
717 Bool_t accept=kTRUE;
718 for(Int_t ilay=0; ilay<6; ilay++){
719 if(fRequirePoint[ilay]){
720 Int_t mask = 1<<ilay;
721 if(!(clumap2 & mask)){
722 accept=kFALSE;
723 break;
724 }
725 }
726 }
727 if(!accept) continue;
728 }
729 // Sort
730 for(Int_t i=0;i<12;i++){
731 for(Int_t j=i+1;j<12;j++){
732 if(idxSA[j]>idxSA[i]){
733 Int_t tmp=idxSA[j];
734 idxSA[j]=idxSA[i];
735 idxSA[i]=tmp;
736 }
737 }
738 }
739 Int_t match=0;
740 for(Int_t icl=0; icl<ncls; icl++){
741 if(idxSA[icl]!=idxMI[icl]){
742 match=0;
743 break;
744 }
745 else match++;
746 }
747 if(match==ncls && match>0){
748 matched=kTRUE;
749 Float_t pt1=track->Pt();
750 Float_t pt2=track2->Pt();
751 Float_t ptm1=track->OneOverPt();
752 Float_t ptm2=track2->OneOverPt();
753 fHistPtResid[hadronSpecie]->Fill(pt1,pt2-pt1);
754 fHistPtRelResid[hadronSpecie]->Fill(pt1,(pt2-pt1)/pt1);
755 fHistInvPtResid[hadronSpecie]->Fill(pt1,ptm2-ptm1);
2afdea69 756 fHistInvPtRelResid[hadronSpecie]->Fill(pt1,(ptm2-ptm1)/ptm1);
757 Float_t deltaphi=track->Phi()-track2->Phi();
758 if(deltaphi<-TMath::Pi()) deltaphi+=2*TMath::Pi();
759 if(deltaphi>TMath::Pi()) deltaphi-=2*TMath::Pi();
760 fHistPhiResid->Fill(pt1,deltaphi);
045571a6 761 }
762 }
763 }
764
765 PostData(1,fOutput);
766
767}
768//______________________________________________________________________________
769void AliAnalysisTaskITSsaTracks::Terminate(Option_t */*option*/)
770{
771 // Terminate analysis
772 fOutput = dynamic_cast<TList*> (GetOutputData(1));
773 if (!fOutput) {
774 printf("ERROR: fOutput not available\n");
775 return;
776 }
2afdea69 777 fNtupleTracks = dynamic_cast<TNtuple*>(fOutput->FindObject("ntupleTracks"));
045571a6 778 fHistNEvents= dynamic_cast<TH1F*>(fOutput->FindObject("hNEvents"));
045571a6 779
2afdea69 780 TString spname[3]={"Pion","Kaon","Proton"};
045571a6 781
782 for(Int_t iSpec=0; iSpec<kNspecies; iSpec++){
783 fHistPtTPCITS[iSpec]= dynamic_cast<TH1F*>(fOutput->FindObject(Form("hPtTPCITS%s",spname[iSpec].Data())));
784 fHistPtITSsa[iSpec]= dynamic_cast<TH1F*>(fOutput->FindObject(Form("hPtITSsa%s",spname[iSpec].Data())));
785 fHistPtITSpureSA[iSpec]= dynamic_cast<TH1F*>(fOutput->FindObject(Form("hPtITSpureSA%s",spname[iSpec].Data())));
786
787 fHistEtaPhiTPCITS[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hEtaPhiTPCITS%s",spname[iSpec].Data())));
788 fHistEtaPhiITSsa[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hEtaPhiITSsa%s",spname[iSpec].Data())));
789 fHistEtaPhiITSpureSA[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hEtaPhiITSpureSA%s",spname[iSpec].Data())));
790
791 fHistNcluTPCITS[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hNcluTPCITS%s",spname[iSpec].Data())));
792 fHistNcluITSsa[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hNcluITSsa%s",spname[iSpec].Data())));
793 fHistNcluITSpureSA[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hNcluITSpureSA%s",spname[iSpec].Data())));
794
795 fHistd0rphiITSpureSA[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hd0rphiITSpureSA%s",spname[iSpec].Data())));
796 fHistd0zITSpureSA[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hd0zITSpureSA%s",spname[iSpec].Data())));
797
798 fHistCluInLayTPCITS[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hCluInLayTPCITS%s",spname[iSpec].Data())));
799 fHistCluInLayITSsa[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hCluInLayITSsa%s",spname[iSpec].Data())));
800 fHistCluInLayITSpureSA[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hCluInLayITSpureSA%s",spname[iSpec].Data())));
801
802 fHistOuterLayITSpureSA[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hOuterLayITSpureSA%s",spname[iSpec].Data())));
803
804 fHistPtResid[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hPtResid%s",spname[iSpec].Data())));
805 fHistPtRelResid[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hPtRelResid%s",spname[iSpec].Data())));
806 fHistInvPtResid[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hInvPtResid%s",spname[iSpec].Data())));
807 fHistInvPtRelResid[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hInvPtRelResid%s",spname[iSpec].Data())));
808 fHistMCPtResid[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hMCPtResid%s",spname[iSpec].Data())));
809 fHistMCPtRelResid[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hMCPtRelResid%s",spname[iSpec].Data())));
810 fHistMCInvPtResid[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hMCInvPtResid%s",spname[iSpec].Data())));
811 fHistMCInvPtRelResid[iSpec]= dynamic_cast<TH2F*>(fOutput->FindObject(Form("hMCInvPtRelResid%s",spname[iSpec].Data())));
812
813 }
2afdea69 814 fHistMCPhiResid= dynamic_cast<TH2F*>(fOutput->FindObject("hMCPhiResid"));
815 fHistPhiResid= dynamic_cast<TH2F*>(fOutput->FindObject("hPhiResid"));
045571a6 816 return;
817}
818
819
820
821
822