]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - PWGJE/AliPWG4HighPtTrackQA.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
Extend MC truth information with z-Vertex and pileup
[u/mrichter/AliRoot.git] / PWGJE / AliPWG4HighPtTrackQA.cxx
... / ...
CommitLineData
1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16//-----------------------------------------------------------------------
17// This class stores QA variables as function of pT for different type
18// of tracks and track selection criteria
19// Output: Histograms for different set of cuts
20//-----------------------------------------------------------------------
21// Author : M. Verweij - UU
22//-----------------------------------------------------------------------
23
24#ifndef ALIPWG4HIGHPTTRACKQA_CXX
25#define ALIPWG4HIGHPTTRACKQA_CXX
26
27#include "AliPWG4HighPtTrackQA.h"
28
29#include "TH1.h"
30#include "TH2.h"
31#include "TH3.h"
32#include "TProfile.h"
33#include "TList.h"
34#include "TFile.h"
35#include "TChain.h"
36#include "TH3F.h"
37#include "TKey.h"
38#include "TSystem.h"
39#include "TBits.h"
40
41#include "AliAnalysisManager.h"
42#include "AliESDInputHandler.h"
43#include "AliMCEvent.h"
44#include "AliMCEventHandler.h"
45#include "AliStack.h"
46#include "AliESDtrack.h"
47#include "AliESDtrackCuts.h"
48#include "AliExternalTrackParam.h"
49#include "AliLog.h"
50#include "AliGenPythiaEventHeader.h"
51#include "AliGenCocktailEventHeader.h"
52#include "AliCentrality.h"
53#include "AliAODVertex.h"
54#include "AliAODEvent.h"
55
56using namespace std; //required for resolving the 'cout' symbol
57
58ClassImp(AliPWG4HighPtTrackQA)
59
60AliPWG4HighPtTrackQA::AliPWG4HighPtTrackQA()
61: AliAnalysisTaskSE(),
62 fDataType(kESD),
63 fEvent(0x0),
64 fESD(0x0),
65 fVtx(0x0),
66 fVtxAOD(0x0),
67 fTrackCuts(0x0),
68 fTrackCutsITSLoose(0x0),
69 fTrackCutsTPConly(0x0),
70 fTrackType(0),
71 fFilterMask(0),
72 fIncludeNoITS(kFALSE),
73 fSigmaConstrainedMax(-1.),
74 fPtMax(100.),
75 fIsPbPb(0),
76 fCentClass(10),
77 fNVariables(26),
78 fVariables(0x0),
79 fITSClusterMap(0),
80 fAvgTrials(1),
81 fNEventAll(0),
82 fNEventSel(0),
83 fNEventReject(0),
84 fh1Centrality(0x0),
85 fh1Xsec(0),
86 fh1Trials(0),
87 fh1PtHard(0),
88 fh1PtHardTrials(0),
89 fh1NTracksAll(0x0),
90 fh1NTracksReject(0x0),
91 fh1NTracksSel(0x0),
92 fPtAll(0),
93 fPtSel(0),
94 fPtPhi(0x0),
95 fPtEta(0x0),
96 fPtEtaPhi(0x0),
97 fPtDCA2D(0x0),
98 fPtDCAZ(0x0),
99 fPtNClustersTPC(0x0),
100 fPtNClustersTPCPhi(0x0),
101 fPtNClustersTPCIter1(0x0),
102 fPtNClustersTPCIter1Phi(0x0),
103 fPtNClustersTPCShared(0x0),
104 fPtNClustersTPCSharedFrac(0x0),
105 fPtNPointITS(0x0),
106 fPtNPointITSPhi(0x0),
107 fPtChi2C(0x0),
108 fPtNSigmaToVertex(0x0),
109 fPtRelUncertainty1Pt(0x0),
110 fPtRelUncertainty1PtNClus(0x0),
111 fPtRelUncertainty1PtNClusIter1(0x0),
112 fPtRelUncertainty1PtNPointITS(0x0),
113 fPtRelUncertainty1PtITSClusterMap(0x0),
114 fPtRelUncertainty1PtChi2(0x0),
115 fPtRelUncertainty1PtChi2Iter1(0x0),
116 fPtRelUncertainty1PtPhi(0x0),
117 fPtChi2PerClusterTPC(0x0),
118 fPtChi2PerClusterTPCIter1(0x0),
119 fPtNCrossedRows(0x0),
120 fPtNCrossedRowsPhi(0x0),
121 fPtNCrossedRowsNClusFPhi(0x0),
122 fPtNCrRNCrRNClusF(0x0),
123 fPtNCrossedRowsFit(0x0),
124 fPtNCrossedRowsFitPhi(0x0),
125 fPtNCrossedRowsNClusFFitPhi(0x0),
126 fNCrossedRowsNCrossedRowsFit(0x0),
127 fNClustersNCrossedRows(0x0),
128 fNClustersNCrossedRowsFit(0x0),
129 fPtNClustersNClustersFitMap(0x0),
130 fPtRelUncertainty1PtNCrossedRows(0x0),
131 fPtRelUncertainty1PtNCrossedRowsFit(0x0),
132 fPtChi2Gold(0x0),
133 fPtChi2GGC(0x0),
134 fPtChi2GoldPhi(0x0),
135 fPtChi2GGCPhi(0x0),
136 fChi2GoldChi2GGC(0x0),
137 fPtChi2ITSPhi(0x0),
138 fPtSigmaY2(0x0),
139 fPtSigmaZ2(0x0),
140 fPtSigmaSnp2(0x0),
141 fPtSigmaTgl2(0x0),
142 fPtSigma1Pt2(0x0),
143 fProfPtSigmaY2(0x0),
144 fProfPtSigmaZ2(0x0),
145 fProfPtSigmaSnp2(0x0),
146 fProfPtSigmaTgl2(0x0),
147 fProfPtSigma1Pt2(0x0),
148 fProfPtSigma1Pt(0x0),
149 fProfPtPtSigma1Pt(0x0),
150 fHistList(0)
151{
152 //
153 // Constructor
154 //
155 SetNVariables(26);
156
157 fPtBinEdges[0][0] = 10.;
158 fPtBinEdges[0][1] = 1.;
159 fPtBinEdges[1][0] = 20.;
160 fPtBinEdges[1][1] = 2.;
161 fPtBinEdges[2][0] = 100.;
162 fPtBinEdges[2][1] = 5.;
163
164}
165//________________________________________________________________________
166AliPWG4HighPtTrackQA::AliPWG4HighPtTrackQA(const char *name):
167 AliAnalysisTaskSE(name),
168 fDataType(kESD),
169 fEvent(0x0),
170 fESD(0x0),
171 fVtx(0x0),
172 fVtxAOD(0x0),
173 fTrackCuts(0x0),
174 fTrackCutsITSLoose(0x0),
175 fTrackCutsTPConly(0x0),
176 fTrackType(0),
177 fFilterMask(0),
178 fIncludeNoITS(kFALSE),
179 fSigmaConstrainedMax(-1.),
180 fPtMax(100.),
181 fIsPbPb(0),
182 fCentClass(10),
183 fNVariables(26),
184 fVariables(0x0),
185 fITSClusterMap(0),
186 fAvgTrials(1),
187 fNEventAll(0),
188 fNEventSel(0),
189 fNEventReject(0),
190 fh1Centrality(0x0),
191 fh1Xsec(0),
192 fh1Trials(0),
193 fh1PtHard(0),
194 fh1PtHardTrials(0),
195 fh1NTracksAll(0x0),
196 fh1NTracksReject(0x0),
197 fh1NTracksSel(0x0),
198 fPtAll(0),
199 fPtSel(0),
200 fPtPhi(0x0),
201 fPtEta(0x0),
202 fPtEtaPhi(0x0),
203 fPtDCA2D(0x0),
204 fPtDCAZ(0x0),
205 fPtNClustersTPC(0x0),
206 fPtNClustersTPCPhi(0x0),
207 fPtNClustersTPCIter1(0x0),
208 fPtNClustersTPCIter1Phi(0x0),
209 fPtNClustersTPCShared(0x0),
210 fPtNClustersTPCSharedFrac(0x0),
211 fPtNPointITS(0x0),
212 fPtNPointITSPhi(0x0),
213 fPtChi2C(0x0),
214 fPtNSigmaToVertex(0x0),
215 fPtRelUncertainty1Pt(0x0),
216 fPtRelUncertainty1PtNClus(0x0),
217 fPtRelUncertainty1PtNClusIter1(0x0),
218 fPtRelUncertainty1PtNPointITS(0x0),
219 fPtRelUncertainty1PtITSClusterMap(0x0),
220 fPtRelUncertainty1PtChi2(0x0),
221 fPtRelUncertainty1PtChi2Iter1(0x0),
222 fPtRelUncertainty1PtPhi(0x0),
223 fPtChi2PerClusterTPC(0x0),
224 fPtChi2PerClusterTPCIter1(0x0),
225 fPtNCrossedRows(0x0),
226 fPtNCrossedRowsPhi(0x0),
227 fPtNCrossedRowsNClusFPhi(0x0),
228 fPtNCrRNCrRNClusF(0x0),
229 fPtNCrossedRowsFit(0x0),
230 fPtNCrossedRowsFitPhi(0x0),
231 fPtNCrossedRowsNClusFFitPhi(0x0),
232 fNCrossedRowsNCrossedRowsFit(0x0),
233 fNClustersNCrossedRows(0x0),
234 fNClustersNCrossedRowsFit(0x0),
235 fPtNClustersNClustersFitMap(0x0),
236 fPtRelUncertainty1PtNCrossedRows(0x0),
237 fPtRelUncertainty1PtNCrossedRowsFit(0x0),
238 fPtChi2Gold(0x0),
239 fPtChi2GGC(0x0),
240 fPtChi2GoldPhi(0x0),
241 fPtChi2GGCPhi(0x0),
242 fChi2GoldChi2GGC(0x0),
243 fPtChi2ITSPhi(0x0),
244 fPtSigmaY2(0x0),
245 fPtSigmaZ2(0x0),
246 fPtSigmaSnp2(0x0),
247 fPtSigmaTgl2(0x0),
248 fPtSigma1Pt2(0x0),
249 fProfPtSigmaY2(0x0),
250 fProfPtSigmaZ2(0x0),
251 fProfPtSigmaSnp2(0x0),
252 fProfPtSigmaTgl2(0x0),
253 fProfPtSigma1Pt2(0x0),
254 fProfPtSigma1Pt(0x0),
255 fProfPtPtSigma1Pt(0x0),
256 fHistList(0)
257{
258 //
259 // Constructor. Initialization of Inputs and Outputs
260 //
261 AliDebug(2,Form("AliPWG4HighPtTrackQA Calling Constructor"));
262
263 SetNVariables(26);
264
265 fPtBinEdges[0][0] = 10.;
266 fPtBinEdges[0][1] = 1.;
267 fPtBinEdges[1][0] = 20.;
268 fPtBinEdges[1][1] = 2.;
269 fPtBinEdges[2][0] = 100.;
270 fPtBinEdges[2][1] = 5.;
271
272 // Input slot #0 works with a TChain ESD
273 DefineInput(0, TChain::Class());
274 // Output slot #1 write into a TList
275 DefineOutput(1, TList::Class());
276}
277
278//________________________________________________________________________
279void AliPWG4HighPtTrackQA::SetPtBinEdges(Int_t region, Double_t ptmax, Double_t ptBinWidth)
280{
281 //
282 // Set variable bin sizes for pT axis in histos
283 //
284
285 if(region<3) {
286 fPtBinEdges[region][0] = ptmax;
287 fPtBinEdges[region][1] = ptBinWidth;
288 }
289 else {
290 AliError("Only 3 regions alowed. Use region 0/1/2\n");
291 return;
292 }
293
294}
295
296//________________________________________________________________________
297void AliPWG4HighPtTrackQA::UserCreateOutputObjects()
298{
299 //Create output objects
300 AliDebug(2,Form(">> AliPWG4HighPtTrackQA::UserCreateOutputObjects \n"));
301
302 Bool_t oldStatus = TH1::AddDirectoryStatus();
303 TH1::AddDirectory(kFALSE);
304
305 OpenFile(1);
306 fHistList = new TList();
307 fHistList->SetOwner(kTRUE);
308
309 Float_t fgkPtMin = 0.;
310 // Float_t fgkPtMax = fPtMax;
311
312 //fPtBinEdges[region][0] = ptmax of region ; fPtBinEdges[region][1] = binWidth of region
313 const Float_t ptmin1 = fgkPtMin;
314 const Float_t ptmax1 = fPtBinEdges[0][0];
315 const Float_t ptmin2 = ptmax1 ;
316 const Float_t ptmax2 = fPtBinEdges[1][0];
317 const Float_t ptmin3 = ptmax2 ;
318 const Float_t ptmax3 = fPtBinEdges[2][0];//fgkPtMax;
319 const Int_t nbin11 = (int)((ptmax1-ptmin1)/fPtBinEdges[0][1]);
320 const Int_t nbin12 = (int)((ptmax2-ptmin2)/fPtBinEdges[1][1])+nbin11;
321 const Int_t nbin13 = (int)((ptmax3-ptmin3)/fPtBinEdges[2][1])+nbin12;
322 Int_t fgkNPtBins=nbin13;
323 //Create array with low edges of each bin
324 Double_t *binsPt=new Double_t[fgkNPtBins+1];
325 for(Int_t i=0; i<=fgkNPtBins; i++) {
326 if(i<=nbin11) binsPt[i]=(Double_t)ptmin1 + (ptmax1-ptmin1)/nbin11*(Double_t)i ;
327 if(i<=nbin12 && i>nbin11) binsPt[i]=(Double_t)ptmin2 + (ptmax2-ptmin2)/(nbin12-nbin11)*((Double_t)i-(Double_t)nbin11) ;
328 if(i<=nbin13 && i>nbin12) binsPt[i]=(Double_t)ptmin3 + (ptmax3-ptmin3)/(nbin13-nbin12)*((Double_t)i-(Double_t)nbin12) ;
329 }
330
331 Int_t fgkNPhiBins = 18*6;
332 Float_t kMinPhi = 0.;
333 Float_t kMaxPhi = 2.*TMath::Pi();
334 Double_t *binsPhi = new Double_t[fgkNPhiBins+1];
335 for(Int_t i=0; i<=fgkNPhiBins; i++) binsPhi[i]=(Double_t)kMinPhi + (kMaxPhi-kMinPhi)/fgkNPhiBins*(Double_t)i ;
336
337 Int_t fgkNEtaBins=20;
338 Float_t fgkEtaMin = -1.;
339 Float_t fgkEtaMax = 1.;
340 Double_t *binsEta=new Double_t[fgkNEtaBins+1];
341 for(Int_t i=0; i<=fgkNEtaBins; i++) binsEta[i]=(Double_t)fgkEtaMin + (fgkEtaMax-fgkEtaMin)/fgkNEtaBins*(Double_t)i ;
342
343 Int_t fgkNNClustersTPCBins=80;
344 Float_t fgkNClustersTPCMin = 0.5;
345 Float_t fgkNClustersTPCMax = 160.5;
346 Double_t *binsNClustersTPC=new Double_t[fgkNNClustersTPCBins+1];
347 for(Int_t i=0; i<=fgkNNClustersTPCBins; i++) binsNClustersTPC[i]=(Double_t)fgkNClustersTPCMin + (fgkNClustersTPCMax-fgkNClustersTPCMin)/fgkNNClustersTPCBins*(Double_t)i ;
348
349 Int_t fgkNDCA2DBins=80;
350 Float_t fgkDCA2DMin = -0.2;
351 Float_t fgkDCA2DMax = 0.2;
352 if(fTrackType==1 || fTrackType==2 || fTrackType==4 || fTrackType==7) {
353 fgkDCA2DMin = -2.;
354 fgkDCA2DMax = 2.;
355 }
356 Double_t *binsDCA2D=new Double_t[fgkNDCA2DBins+1];
357 for(Int_t i=0; i<=fgkNDCA2DBins; i++) binsDCA2D[i]=(Double_t)fgkDCA2DMin + (fgkDCA2DMax-fgkDCA2DMin)/fgkNDCA2DBins*(Double_t)i ;
358
359 Int_t fgkNDCAZBins=80;
360 Float_t fgkDCAZMin = -2.;
361 Float_t fgkDCAZMax = 2.;
362 if(fTrackType==1 || fTrackType==2 || fTrackType==4) {
363 fgkDCAZMin = -5.;
364 fgkDCAZMax = 5.;
365 }
366 Double_t *binsDCAZ=new Double_t[fgkNDCAZBins+1];
367 for(Int_t i=0; i<=fgkNDCAZBins; i++) binsDCAZ[i]=(Double_t)fgkDCAZMin + (fgkDCAZMax-fgkDCAZMin)/fgkNDCAZBins*(Double_t)i ;
368
369 Int_t fgkNNPointITSBins=9;
370 Float_t fgkNPointITSMin = -0.5;
371 Float_t fgkNPointITSMax = 8.5;
372 Double_t *binsNPointITS=new Double_t[fgkNNPointITSBins+1];
373 for(Int_t i=0; i<=fgkNNPointITSBins; i++) binsNPointITS[i]=(Double_t)fgkNPointITSMin + (fgkNPointITSMax-fgkNPointITSMin)/fgkNNPointITSBins*(Double_t)i ;
374
375 Int_t fgkNITSClusterMapBins=65;
376 Float_t fgkITSClusterMapMin = -0.5;
377 Float_t fgkITSClusterMapMax = 64.5;
378 Double_t *binsITSClusterMap=new Double_t[fgkNITSClusterMapBins+1];
379 for(Int_t i=0; i<=fgkNITSClusterMapBins; i++) binsITSClusterMap[i]=(Double_t)fgkITSClusterMapMin + (fgkITSClusterMapMax-fgkITSClusterMapMin)/fgkNITSClusterMapBins*(Double_t)i ;
380
381
382 Int_t fgkNNSigmaToVertexBins=9;
383 Float_t fgkNSigmaToVertexMin = 0.;
384 Float_t fgkNSigmaToVertexMax = 9.;
385 Double_t *binsNSigmaToVertex=new Double_t[fgkNNSigmaToVertexBins+1];
386 for(Int_t i=0; i<=fgkNNSigmaToVertexBins; i++) binsNSigmaToVertex[i]=(Double_t)fgkNSigmaToVertexMin + (fgkNSigmaToVertexMax-fgkNSigmaToVertexMin)/fgkNNSigmaToVertexBins*(Double_t)i ;
387
388 Int_t fgkNChi2CBins=10;
389 // Float_t fgkChi2CMin = 0.;
390 // Float_t fgkChi2CMax = 100.; //10 sigma
391 Double_t *binsChi2C=new Double_t[fgkNChi2CBins+1];
392 for(Int_t i=0; i<=fgkNChi2CBins; i++) binsChi2C[i] = (Double_t)i * (Double_t)i;
393
394 Float_t fgkRel1PtUncertaintyMin = 0.;
395 Float_t fgkRel1PtUncertaintyMax = 1.;
396 Float_t binEdgeRel1PtUncertainty1= 0.3;
397 Int_t fgkNRel1PtUncertaintyBins1 = 45;
398 Float_t binWidthRel1PtUncertainty1 = (binEdgeRel1PtUncertainty1-fgkRel1PtUncertaintyMin)/((Float_t)fgkNRel1PtUncertaintyBins1);
399 Int_t fgkNRel1PtUncertaintyBins2 = 35;
400 Float_t binWidthRel1PtUncertainty2 = (fgkRel1PtUncertaintyMax-binEdgeRel1PtUncertainty1)/((Float_t)fgkNRel1PtUncertaintyBins2);
401 Int_t fgkNRel1PtUncertaintyBins = fgkNRel1PtUncertaintyBins1 + fgkNRel1PtUncertaintyBins2;
402
403 Double_t *binsRel1PtUncertainty=new Double_t[fgkNRel1PtUncertaintyBins+1];
404 for(Int_t i=0; i<=fgkNRel1PtUncertaintyBins; i++) {
405 if(i<=fgkNRel1PtUncertaintyBins1)
406 binsRel1PtUncertainty[i]=(Double_t)fgkRel1PtUncertaintyMin + (Double_t)binWidthRel1PtUncertainty1*(Double_t)i ;
407 if(i<=fgkNRel1PtUncertaintyBins && i>fgkNRel1PtUncertaintyBins1)
408 binsRel1PtUncertainty[i]=(Double_t)binEdgeRel1PtUncertainty1 + (Double_t)binWidthRel1PtUncertainty2*(Double_t)(i-fgkNRel1PtUncertaintyBins1);
409 }
410
411 Int_t fgkNUncertainty1PtBins = 30;
412 Float_t fgkUncertainty1PtMin = 0.;
413 Float_t fgkUncertainty1PtMax = 0.1;
414 if(fTrackType==1 || fTrackType==2 || fTrackType==4)
415 fgkUncertainty1PtMax = 0.2;
416 Double_t *binsUncertainty1Pt=new Double_t[fgkNUncertainty1PtBins+1];
417 for(Int_t i=0; i<=fgkNUncertainty1PtBins; i++) binsUncertainty1Pt[i]=(Double_t)fgkUncertainty1PtMin + (fgkUncertainty1PtMax-fgkUncertainty1PtMin)/fgkNUncertainty1PtBins*(Double_t)i ;
418
419 Float_t fgkChi2PerClusMin = 0.;
420 Float_t fgkChi2PerClusMax = 4.;
421 Int_t fgkNChi2PerClusBins = (int)(fgkChi2PerClusMax*10.);
422 Double_t *binsChi2PerClus=new Double_t[fgkNChi2PerClusBins+1];
423 for(Int_t i=0; i<=fgkNChi2PerClusBins; i++) binsChi2PerClus[i]=(Double_t)fgkChi2PerClusMin + (fgkChi2PerClusMax-fgkChi2PerClusMin)/fgkNChi2PerClusBins*(Double_t)i ;
424
425 Int_t fgkNCrossedRowsNClusFBins = 45;
426 Float_t fgkNCrossedRowsNClusFMin = 0.;
427 Float_t fgkNCrossedRowsNClusFMax = 1.5;
428 Double_t *binsNCrossedRowsNClusF=new Double_t[fgkNCrossedRowsNClusFBins+1];
429 for(Int_t i=0; i<=fgkNCrossedRowsNClusFBins; i++) binsNCrossedRowsNClusF[i]=(Double_t)fgkNCrossedRowsNClusFMin + (fgkNCrossedRowsNClusFMax-fgkNCrossedRowsNClusFMin)/fgkNCrossedRowsNClusFBins*(Double_t)i ;
430
431 Float_t fgk1PtMin = 0.;
432 Float_t fgk1PtMax = 6.;
433 Float_t binEdge1Pt1 = 1.;
434 Float_t binWidth1Pt1 = 0.05;
435 Int_t fgkN1PtBins1 = (int)((binEdge1Pt1-fgk1PtMin)/binWidth1Pt1);
436 Float_t binWidth1Pt2 = 0.1;
437 Int_t fgkN1PtBins2 = (int)((fgk1PtMax-binEdge1Pt1)/binWidth1Pt2);
438 Int_t fgkN1PtBins = fgkN1PtBins1+fgkN1PtBins2;
439 Double_t *bins1Pt=new Double_t[fgkN1PtBins+1];
440
441 for(Int_t i=0; i<=fgkN1PtBins; i++) {
442 if(i<=fgkN1PtBins1)
443 bins1Pt[i]=(Double_t)fgk1PtMin + (Double_t)(binEdge1Pt1-fgk1PtMin)/(Double_t)fgkN1PtBins1*(Double_t)i;
444 if(i<=fgkN1PtBins && i>fgkN1PtBins1)
445 bins1Pt[i]=(Double_t)binEdge1Pt1 + (Double_t)(fgk1PtMax-binEdge1Pt1)/(Double_t)fgkN1PtBins2*(Double_t)(i-fgkN1PtBins1);
446 }
447
448 Int_t fgkNSigmaY2Bins = 50;
449 Float_t fgkSigmaY2Min = 0.;
450 Float_t fgkSigmaY2Max = 1.;
451 if(fTrackType==1) fgkSigmaY2Max = 4.;
452 if(fTrackType==2 || fTrackType==4) fgkSigmaY2Max = 0.1;
453 Double_t *binsSigmaY2=new Double_t[fgkNSigmaY2Bins+1];
454 for(Int_t i=0; i<=fgkNSigmaY2Bins; i++) binsSigmaY2[i]=(Double_t)fgkSigmaY2Min + (fgkSigmaY2Max-fgkSigmaY2Min)/fgkNSigmaY2Bins*(Double_t)i ;
455
456 Int_t fgkNSigmaZ2Bins = 50;
457 Float_t fgkSigmaZ2Min = 0.;
458 Float_t fgkSigmaZ2Max = 0.4;
459 Double_t *binsSigmaZ2=new Double_t[fgkNSigmaZ2Bins+1];
460 for(Int_t i=0; i<=fgkNSigmaZ2Bins; i++) binsSigmaZ2[i]=(Double_t)fgkSigmaZ2Min + (fgkSigmaZ2Max-fgkSigmaZ2Min)/fgkNSigmaZ2Bins*(Double_t)i ;
461
462 Int_t fgkNSigmaSnp2Bins = 50;
463 Float_t fgkSigmaSnp2Min = 0.;
464 Float_t fgkSigmaSnp2Max = 0.05;
465 if(fTrackType==1) fgkSigmaSnp2Max = 0.2;
466 if(fTrackType==2 || fTrackType==4) fgkSigmaSnp2Max = 0.1;
467 Double_t *binsSigmaSnp2=new Double_t[fgkNSigmaSnp2Bins+1];
468 for(Int_t i=0; i<=fgkNSigmaSnp2Bins; i++) binsSigmaSnp2[i]=(Double_t)fgkSigmaSnp2Min + (fgkSigmaSnp2Max-fgkSigmaSnp2Min)/fgkNSigmaSnp2Bins*(Double_t)i ;
469
470 Int_t fgkNSigmaTgl2Bins = 50;
471 Float_t fgkSigmaTgl2Min = 0.;
472 Float_t fgkSigmaTgl2Max = 0.1;
473 if(fTrackType==1) fgkSigmaTgl2Max = 0.2;
474 if(fTrackType==2 || fTrackType==4) fgkSigmaTgl2Max = 0.1;
475 Double_t *binsSigmaTgl2=new Double_t[fgkNSigmaTgl2Bins+1];
476 for(Int_t i=0; i<=fgkNSigmaTgl2Bins; i++) binsSigmaTgl2[i]=(Double_t)fgkSigmaTgl2Min + (fgkSigmaTgl2Max-fgkSigmaTgl2Min)/fgkNSigmaTgl2Bins*(Double_t)i ;
477
478 Int_t fgkNSigma1Pt2Bins = 50;
479 Float_t fgkSigma1Pt2Min = 0.;
480 Float_t fgkSigma1Pt2Max = 1.;
481 Double_t *binsSigma1Pt2=new Double_t[fgkNSigma1Pt2Bins+1];
482 for(Int_t i=0; i<=fgkNSigma1Pt2Bins; i++) binsSigma1Pt2[i]=(Double_t)fgkSigma1Pt2Min + (fgkSigma1Pt2Max-fgkSigma1Pt2Min)/fgkNSigma1Pt2Bins*(Double_t)i ;
483
484
485 fNEventAll = new TH1F("fNEventAll","NEventAll",1,-0.5,0.5);
486 fHistList->Add(fNEventAll);
487 fNEventSel = new TH1F("fNEventSel","NEvent Selected for analysis",1,-0.5,0.5);
488 fHistList->Add(fNEventSel);
489 fNEventReject = new TH1F("fNEventReject","Reason events are rejectected for analysis",20,0,20);
490 //Set labels
491 fNEventReject->Fill("noESD",0);
492 fNEventReject->Fill("Trigger",0);
493 fNEventReject->Fill("NTracks<2",0);
494 fNEventReject->Fill("noVTX",0);
495 fNEventReject->Fill("VtxStatus",0);
496 fNEventReject->Fill("NCont<2",0);
497 fNEventReject->Fill("ZVTX>10",0);
498 fNEventReject->Fill("cent",0);
499 fNEventReject->Fill("cent>90",0);
500 fHistList->Add(fNEventReject);
501
502 fh1Centrality = new TH1F("fh1Centrality","fh1Centrality; Centrality %",100,0,100);
503 fHistList->Add(fh1Centrality);
504
505 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
506 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
507 fHistList->Add(fh1Xsec);
508
509 fh1Trials = new TH1F("fh1Trials","trials root file",1,0,1);
510 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
511 fHistList->Add(fh1Trials);
512
513 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
514 fHistList->Add(fh1PtHard);
515 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
516 fHistList->Add(fh1PtHardTrials);
517
518 fh1NTracksAll = new TH1F("fh1NTracksAll","fh1NTracksAll",1,-0.5,0.5);
519 fHistList->Add(fh1NTracksAll);
520
521 fh1NTracksReject = new TH1F("fh1NTracksReject","fh1NTracksReject",1,-0.5,0.5);
522 fh1NTracksReject->Fill("noHybridTrack",0);
523 fh1NTracksReject->Fill("noITSrefit",0);
524 fh1NTracksReject->Fill("noESDtrack",0);
525 fh1NTracksReject->Fill("noTPCInner",0);
526 fh1NTracksReject->Fill("FillTPC",0);
527 fh1NTracksReject->Fill("noTPConly",0);
528 fh1NTracksReject->Fill("relate",0);
529 fh1NTracksReject->Fill("trackCuts",0);
530 fh1NTracksReject->Fill("laser",0);
531 fh1NTracksReject->Fill("chi2",0);
532 fHistList->Add(fh1NTracksReject);
533
534 fh1NTracksSel = new TH1F("fh1NTracksSel","fh1NTracksSel",1,-0.5,0.5);
535 fHistList->Add(fh1NTracksSel);
536
537 fPtAll = new TH1F("fPtAll","PtAll",fgkNPtBins, binsPt);
538 fHistList->Add(fPtAll);
539 fPtSel = new TH1F("fPtSel","PtSel",fgkNPtBins, binsPt);
540 fHistList->Add(fPtSel);
541
542 fPtPhi = new TH2F("fPtPhi","fPtPhi",fgkNPtBins,binsPt,fgkNPhiBins,binsPhi);
543 fHistList->Add(fPtPhi);
544
545 fPtEta = new TH2F("fPtEta","fPtEta",fgkNPtBins,binsPt,fgkNEtaBins,binsEta);
546 fHistList->Add(fPtEta);
547
548 fPtEtaPhi = new TH3F("fPtEtaPhi","fPtEtaPhi",fgkNPtBins,binsPt,fgkNEtaBins,binsEta,fgkNPhiBins,binsPhi);
549 fHistList->Add(fPtEtaPhi);
550
551 fPtDCA2D = new TH2F("fPtDCA2D","fPtDCA2D",fgkNPtBins,binsPt,fgkNDCA2DBins,binsDCA2D);
552 fHistList->Add(fPtDCA2D);
553
554 fPtDCAZ = new TH2F("fPtDCAZ","fPtDCAZ",fgkNPtBins,binsPt,fgkNDCAZBins,binsDCAZ);
555 fHistList->Add(fPtDCAZ);
556
557 fPtNClustersTPC = new TH2F("fPtNClustersTPC","fPtNClustersTPC",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC);
558 fHistList->Add(fPtNClustersTPC);
559
560 fPtNClustersTPCPhi = new TH2F("fPtNClustersTPCPhi","fPtNClustersTPCPhi",fgkNPhiBins,binsPhi,fgkNNClustersTPCBins,binsNClustersTPC);
561 fHistList->Add(fPtNClustersTPCPhi);
562
563 fPtNClustersTPCIter1 = new TH2F("fPtNClustersTPCIter1","fPtNClustersTPCIter1",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC);
564 fHistList->Add(fPtNClustersTPCIter1);
565
566 fPtNClustersTPCIter1Phi = new TH3F("fPtNClustersTPCIter1Phi","fPtNClustersTPCIter1Phi",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC,fgkNPhiBins,binsPhi);
567 fHistList->Add(fPtNClustersTPCIter1Phi);
568
569 fPtNClustersTPCShared = new TH2F("fPtNClustersTPCShared","fPtNClustersTPCShared",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC);
570 fHistList->Add(fPtNClustersTPCShared);
571
572 fPtNClustersTPCSharedFrac = new TH2F("fPtNClustersTPCSharedFrac","fPtNClustersTPCSharedFrac",fgkNPtBins,binsPt,fgkNSigma1Pt2Bins,binsSigma1Pt2);
573 fHistList->Add(fPtNClustersTPCSharedFrac);
574
575 fPtNPointITS = new TH2F("fPtNPointITS","fPtNPointITS",fgkNPtBins,binsPt,fgkNNPointITSBins,binsNPointITS);
576 fHistList->Add(fPtNPointITS);
577
578 fPtNPointITSPhi = new TH3F("fPtNPointITSPhi","fPtNPointITSPhi",fgkNPtBins,binsPt,fgkNNPointITSBins,binsNPointITS,fgkNPhiBins,binsPhi);
579 fHistList->Add(fPtNPointITSPhi);
580
581 fPtChi2C = new TH2F("fPtChi2C","fPtChi2C",fgkNPtBins,binsPt,fgkNChi2CBins,binsChi2C);
582 fHistList->Add(fPtChi2C);
583
584 fPtNSigmaToVertex = new TH2F("fPtNSigmaToVertex","fPtNSigmaToVertex",fgkNPtBins,binsPt,fgkNNSigmaToVertexBins,binsNSigmaToVertex);
585 fHistList->Add(fPtNSigmaToVertex);
586
587 fPtRelUncertainty1Pt = new TH2F("fPtRelUncertainty1Pt","fPtRelUncertainty1Pt",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty);
588 fHistList->Add(fPtRelUncertainty1Pt);
589
590 fPtRelUncertainty1PtNClus = new TH3F("fPtRelUncertainty1PtNClus","fPtRelUncertainty1PtNClus",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNNClustersTPCBins,binsNClustersTPC);
591 fHistList->Add(fPtRelUncertainty1PtNClus);
592
593 fPtRelUncertainty1PtNClusIter1 = new TH3F("fPtRelUncertainty1PtNClusIter1","fPtRelUncertainty1PtNClusIter1",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNNClustersTPCBins,binsNClustersTPC);
594 fHistList->Add(fPtRelUncertainty1PtNClusIter1);
595
596 fPtRelUncertainty1PtNPointITS = new TH3F("fPtRelUncertainty1PtNPointITS","fPtRelUncertainty1PtNPointITS",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNNPointITSBins,binsNPointITS);
597 fHistList->Add(fPtRelUncertainty1PtNPointITS);
598
599 fPtRelUncertainty1PtITSClusterMap = new TH3F("fPtRelUncertainty1PtITSClusterMap","fPtRelUncertainty1PtITSClusterMap",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNITSClusterMapBins,binsITSClusterMap);
600 fHistList->Add(fPtRelUncertainty1PtITSClusterMap);
601
602 fPtRelUncertainty1PtChi2 = new TH3F("fPtRelUncertainty1PtChi2","fPtRelUncertainty1PtChi2",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNChi2PerClusBins,binsChi2PerClus);
603 fHistList->Add(fPtRelUncertainty1PtChi2);
604
605 fPtRelUncertainty1PtChi2Iter1 = new TH3F("fPtRelUncertainty1PtChi2Iter1","fPtRelUncertainty1PtChi2Iter1",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNChi2PerClusBins,binsChi2PerClus);
606 fHistList->Add(fPtRelUncertainty1PtChi2Iter1);
607
608 fPtRelUncertainty1PtPhi = new TH3F("fPtRelUncertainty1PtPhi","fPtRelUncertainty1PtPhi",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNPhiBins,binsPhi);
609 fHistList->Add(fPtRelUncertainty1PtPhi);
610
611 fPtChi2PerClusterTPC = new TH2F("fPtChi2PerClusterTPC","fPtChi2PerClusterTPC",fgkNPtBins,binsPt,fgkNChi2PerClusBins,binsChi2PerClus);
612 fHistList->Add(fPtChi2PerClusterTPC);
613
614 fPtChi2PerClusterTPCIter1 = new TH2F("fPtChi2PerClusterTPCIter1","fPtChi2PerClusterTPCIter1",fgkNPtBins,binsPt,fgkNChi2PerClusBins,binsChi2PerClus);
615 fHistList->Add(fPtChi2PerClusterTPCIter1);
616
617 fPtNCrossedRows = new TH2F("fPtNCrossedRows","fPtNCrossedRows",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC);
618 fHistList->Add(fPtNCrossedRows);
619
620 fPtNCrossedRowsPhi = new TH3F("fPtNCrossedRowsPhi","fPtNCrossedRowsPhi",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC,fgkNPhiBins,binsPhi);
621 fHistList->Add(fPtNCrossedRowsPhi);
622
623 fPtNCrossedRowsNClusFPhi = new TH3F("fPtNCrossedRowsNClusFPhi","fPtNCrossedRowsNClusFPhi",fgkNPtBins,binsPt,fgkNCrossedRowsNClusFBins,binsNCrossedRowsNClusF,fgkNPhiBins,binsPhi);
624 fHistList->Add(fPtNCrossedRowsNClusFPhi);
625
626 fPtNCrRNCrRNClusF = new TH3F("fPtNCrRNCrRNClusF","fPtNCrRNCrRNClusF",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC,fgkNCrossedRowsNClusFBins,binsNCrossedRowsNClusF);
627 fHistList->Add(fPtNCrRNCrRNClusF);
628
629 fPtNCrossedRowsFit = new TH2F("fPtNCrossedRowsFit","fPtNCrossedRowsFit",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC);
630 fHistList->Add(fPtNCrossedRowsFit);
631
632 fPtNCrossedRowsFitPhi = new TH3F("fPtNCrossedRowsFitPhi","fPtNCrossedRowsFitPhi",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC,fgkNPhiBins,binsPhi);
633 fHistList->Add(fPtNCrossedRowsFitPhi);
634
635 fPtNCrossedRowsNClusFFitPhi = new TH3F("fPtNCrossedRowsNClusFFitPhi","fPtNCrossedRowsNClusFFitPhi",fgkNPtBins,binsPt,fgkNCrossedRowsNClusFBins,binsNCrossedRowsNClusF,fgkNPhiBins,binsPhi);
636 fHistList->Add(fPtNCrossedRowsNClusFFitPhi);
637
638 fNCrossedRowsNCrossedRowsFit = new TH2F("fNCrossedRowsNCrossedRowsFit","fNCrossedRowsNCrossedRowsFit",fgkNNClustersTPCBins,binsNClustersTPC,fgkNNClustersTPCBins,binsNClustersTPC);
639 fHistList->Add(fNCrossedRowsNCrossedRowsFit);
640
641 fNClustersNCrossedRows = new TH2F("fNClustersNCrossedRows","fNClustersNCrossedRows",fgkNNClustersTPCBins,binsNClustersTPC,fgkNNClustersTPCBins,binsNClustersTPC);
642 fHistList->Add(fNClustersNCrossedRows);
643
644 fNClustersNCrossedRowsFit = new TH2F("fNClustersNCrossedRowsFit","fNClustersNCrossedRowsFit",fgkNNClustersTPCBins,binsNClustersTPC,fgkNNClustersTPCBins,binsNClustersTPC);
645 fHistList->Add(fNClustersNCrossedRowsFit);
646
647 fPtNClustersNClustersFitMap = new TH3F("fPtNClustersNClustersFitMap","fPtNClustersNClustersFitMap;p_{T};N_{cls};N_{cls}^{fit map}",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC,fgkNNClustersTPCBins,binsNClustersTPC);
648 fHistList->Add(fPtNClustersNClustersFitMap);
649
650 fPtRelUncertainty1PtNCrossedRows = new TH3F("fPtRelUncertainty1PtNCrossedRows","fPtRelUncertainty1PtNCrossedRows",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNNClustersTPCBins,binsNClustersTPC);
651 fHistList->Add(fPtRelUncertainty1PtNCrossedRows);
652
653 fPtRelUncertainty1PtNCrossedRowsFit = new TH3F("fPtRelUncertainty1PtNCrossedRowsFit","fPtRelUncertainty1PtNCrossedRowsFit",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNNClustersTPCBins,binsNClustersTPC);
654 fHistList->Add(fPtRelUncertainty1PtNCrossedRowsFit);
655
656 fPtChi2Gold = new TH2F("fPtChi2Gold","fPtChi2Gold",fgkNPtBins,binsPt,fgkNChi2CBins,binsChi2C);
657 fHistList->Add(fPtChi2Gold);
658
659 fPtChi2GGC = new TH2F("fPtChi2GGC","fPtChi2GGC",fgkNPtBins,binsPt,fgkNChi2CBins,binsChi2C);
660 fHistList->Add(fPtChi2GGC);
661
662 fPtChi2GoldPhi = new TH3F("fPtChi2GoldPhi","fPtChi2GoldPhi",fgkNPtBins,binsPt,fgkNChi2CBins,binsChi2C,fgkNPhiBins,binsPhi);
663 fHistList->Add(fPtChi2GoldPhi);
664
665 fPtChi2GGCPhi = new TH3F("fPtChi2GGCPhi","fPtChi2GGCPhi",fgkNPtBins,binsPt,fgkNChi2CBins,binsChi2C,fgkNPhiBins,binsPhi);
666 fHistList->Add(fPtChi2GGCPhi);
667
668 fChi2GoldChi2GGC = new TH2F("fChi2GoldChi2GGC","fChi2GoldChi2GGC;#chi^{2}_{gold};#chi^{2}_{ggc}",fgkNChi2CBins,binsChi2C,fgkNChi2CBins,binsChi2C);
669 fHistList->Add(fChi2GoldChi2GGC);
670
671 fPtChi2ITSPhi = new TH3F("fPtChi2ITSPhi","fPtChi2ITSPhi;p_{T};#chi^{2}_{ITS};#varphi",fgkNPtBins,binsPt,fgkNChi2CBins,binsChi2C,fgkNPhiBins,binsPhi);
672 fHistList->Add(fPtChi2ITSPhi);
673
674 fPtSigmaY2 = new TH2F("fPtSigmaY2","fPtSigmaY2",fgkN1PtBins,bins1Pt,fgkNSigmaY2Bins,binsSigmaY2);
675 fHistList->Add(fPtSigmaY2);
676
677 fPtSigmaZ2 = new TH2F("fPtSigmaZ2","fPtSigmaZ2",fgkN1PtBins,bins1Pt,fgkNSigmaZ2Bins,binsSigmaZ2);
678 fHistList->Add(fPtSigmaZ2);
679
680 fPtSigmaSnp2 = new TH2F("fPtSigmaSnp2","fPtSigmaSnp2",fgkN1PtBins,bins1Pt,fgkNSigmaSnp2Bins,binsSigmaSnp2);
681 fHistList->Add(fPtSigmaSnp2);
682
683 fPtSigmaTgl2 = new TH2F("fPtSigmaTgl2","fPtSigmaTgl2",fgkN1PtBins,bins1Pt,fgkNSigmaTgl2Bins,binsSigmaTgl2);
684 fHistList->Add(fPtSigmaTgl2);
685
686 fPtSigma1Pt2 = new TH2F("fPtSigma1Pt2","fPtSigma1Pt2",fgkN1PtBins,bins1Pt,fgkNSigma1Pt2Bins,binsSigma1Pt2);
687 fHistList->Add(fPtSigma1Pt2);
688
689 fProfPtSigmaY2 = new TProfile("fProfPtSigmaY2","fProfPtSigmaY2",fgkN1PtBins,bins1Pt);
690 fHistList->Add(fProfPtSigmaY2);
691
692 fProfPtSigmaZ2 = new TProfile("fProfPtSigmaZ2","fProfPtSigmaZ2",fgkN1PtBins,bins1Pt);
693 fHistList->Add(fProfPtSigmaZ2);
694
695 fProfPtSigmaSnp2 = new TProfile("fProfPtSigmaSnp2","fProfPtSigmaSnp2",fgkN1PtBins,bins1Pt);
696 fHistList->Add(fProfPtSigmaSnp2);
697
698 fProfPtSigmaTgl2 = new TProfile("fProfPtSigmaTgl2","fProfPtSigmaTgl2",fgkN1PtBins,bins1Pt);
699 fHistList->Add(fProfPtSigmaTgl2);
700
701 fProfPtSigma1Pt2 = new TProfile("fProfPtSigma1Pt2","fProfPtSigma1Pt2",fgkN1PtBins,bins1Pt);
702 fHistList->Add(fProfPtSigma1Pt2);
703
704 fProfPtSigma1Pt = new TProfile("fProfPtSigma1Pt","fProfPtSigma1Pt;p_{T};#sigma(1/p_{T})",fgkNPtBins,binsPt);
705 fHistList->Add(fProfPtSigma1Pt);
706
707 fProfPtPtSigma1Pt = new TProfile("fProfPtPtSigma1Pt","fProfPtPtSigma1Pt;p_{T};p_{T}#sigma(1/p_{T})",fgkNPtBins,binsPt);
708 fHistList->Add(fProfPtPtSigma1Pt);
709
710 TH1::AddDirectory(oldStatus);
711
712 PostData(1, fHistList);
713
714 if(binsPhi) delete [] binsPhi;
715 if(binsPt) delete [] binsPt;
716 if(binsNClustersTPC) delete [] binsNClustersTPC;
717 if(binsDCA2D) delete [] binsDCA2D;
718 if(binsDCAZ) delete [] binsDCAZ;
719 if(binsNPointITS) delete [] binsNPointITS;
720 if(binsITSClusterMap) delete [] binsITSClusterMap;
721 if(binsNSigmaToVertex) delete [] binsNSigmaToVertex;
722 if(binsChi2C) delete [] binsChi2C;
723 if(binsEta) delete [] binsEta;
724 if(binsRel1PtUncertainty) delete [] binsRel1PtUncertainty;
725 if(binsUncertainty1Pt) delete [] binsUncertainty1Pt;
726 if(binsChi2PerClus) delete [] binsChi2PerClus;
727 if(binsChi2PerClus) delete [] binsNCrossedRowsNClusF;
728 if(bins1Pt) delete [] bins1Pt;
729 if(binsSigmaY2) delete [] binsSigmaY2;
730 if(binsSigmaZ2) delete [] binsSigmaZ2;
731 if(binsSigmaSnp2) delete [] binsSigmaSnp2;
732 if(binsSigmaTgl2) delete [] binsSigmaTgl2;
733 if(binsSigma1Pt2) delete [] binsSigma1Pt2;
734}
735
736//________________________________________________________________________
737Bool_t AliPWG4HighPtTrackQA::SelectEvent()
738{
739 //
740 // Decide if event should be selected for analysis
741 //
742
743 // Checks following requirements:
744 // - fEvent available
745 // - trigger info from AliPhysicsSelection
746 // - MCevent available
747 // - number of reconstructed tracks > 1
748 // - primary vertex reconstructed
749 // - z-vertex < 10 cm
750 // - centrality in case of PbPb
751
752 Bool_t selectEvent = kTRUE;
753
754 //fEvent object available?
755 if (!fEvent) {
756 AliDebug(2,Form("ERROR: fInputEvent not available\n"));
757 fNEventReject->Fill("noAliVEvent",1);
758 selectEvent = kFALSE;
759 return selectEvent;
760 }
761
762 //Check if number of reconstructed tracks is larger than 1
763 if(!fEvent->GetNumberOfTracks() || fEvent->GetNumberOfTracks()<2) {
764 fNEventReject->Fill("NTracks<2",1);
765 selectEvent = kFALSE;
766 return selectEvent;
767 }
768
769 //Check if vertex is reconstructed
770 if(fDataType==kESD&&dynamic_cast<AliESDEvent*>(fEvent)) {
771 fVtx = ((AliESDEvent*)fEvent)->GetPrimaryVertexTracks();
772
773 if (!fVtx || !fVtx->GetStatus())
774 fVtx = ((AliESDEvent*)fEvent)->GetPrimaryVertexSPD();
775
776 if(!fVtx) {
777 fNEventReject->Fill("noVTX",1);
778 selectEvent = kFALSE;
779 return selectEvent;
780 }
781
782 if(!fVtx->GetStatus()) {
783 fNEventReject->Fill("VtxStatus",1);
784 selectEvent = kFALSE;
785 return selectEvent;
786 }
787
788 // Need vertex cut
789 if(fVtx->GetNContributors()<2) {
790 fNEventReject->Fill("NCont<2",1);
791 selectEvent = kFALSE;
792 return selectEvent;
793 }
794
795 //Check if z-vertex < 10 cm
796 double primVtx[3];
797 fVtx->GetXYZ(primVtx);
798 if(TMath::Sqrt(primVtx[0]*primVtx[0] + primVtx[1]*primVtx[1])>1. || TMath::Abs(primVtx[2]>10.)){
799 fNEventReject->Fill("ZVTX>10",1);
800 selectEvent = kFALSE;
801 return selectEvent;
802 }
803 }
804 else if(fDataType==kAOD&&dynamic_cast<AliAODEvent*>(fEvent)) {
805 fVtxAOD = ((AliAODEvent*)fEvent)->GetPrimaryVertex();
806 if(!fVtxAOD) {
807 fNEventReject->Fill("noVTX",1);
808 selectEvent = kFALSE;
809 return selectEvent;
810 }
811
812 // Need vertex cut
813 if(fVtxAOD->GetNContributors()<2) {
814 fNEventReject->Fill("NCont<2",1);
815 selectEvent = kFALSE;
816 return selectEvent;
817 }
818
819 //Check if z-vertex < 10 cm
820 double primVtx[3];
821 fVtxAOD->GetXYZ(primVtx);
822 if(TMath::Sqrt(primVtx[0]*primVtx[0] + primVtx[1]*primVtx[1])>1. || TMath::Abs(primVtx[2]>10.)){
823 fNEventReject->Fill("ZVTX>10",1);
824 selectEvent = kFALSE;
825 return selectEvent;
826 }
827
828 }
829
830 //Centrality selection should only be done in case of PbPb
831 if(IsPbPb()) {
832 Float_t cent = 0.;
833 if(fCentClass!=CalculateCentrality(fEvent) && fCentClass!=10) {
834 fNEventReject->Fill("cent",1);
835 selectEvent = kFALSE;
836 return selectEvent;
837 }
838 else {
839 if(fDataType==kESD) {
840 if(dynamic_cast<AliESDEvent*>(fEvent)->GetCentrality()) {
841 cent = dynamic_cast<AliESDEvent*>(fEvent)->GetCentrality()->GetCentralityPercentile("V0M");
842 }
843 }
844 else if(fDataType==kAOD) {
845 if(dynamic_cast<AliAODEvent*>(fEvent)->GetHeader()->GetCentrality())
846 cent = dynamic_cast<AliAODEvent*>(fEvent)->GetHeader()->GetCentrality();
847 }
848 if(cent>90.) {
849 fNEventReject->Fill("cent>90",1);
850 selectEvent = kFALSE;
851 return selectEvent;
852 }
853 fh1Centrality->Fill(cent);
854 }
855 }
856
857 return selectEvent;
858
859}
860
861//________________________________________________________________________
862Int_t AliPWG4HighPtTrackQA::CalculateCentrality(AliVEvent *ev)
863{
864 //
865 // Get centrality from ESD or AOD
866 //
867
868 if(fDataType==kESD)
869 return CalculateCentrality(dynamic_cast<AliESDEvent*>(ev));
870 else if(fDataType==kAOD)
871 return CalculateCentrality(dynamic_cast<AliAODEvent*>(ev));
872 else
873 return 5;
874}
875
876//________________________________________________________________________
877Int_t AliPWG4HighPtTrackQA::CalculateCentrality(AliESDEvent *esd)
878{
879 //
880 // Get centrality from ESD
881 //
882
883 Float_t cent = -1;
884
885 if(esd){
886 if(esd->GetCentrality()){
887 cent = esd->GetCentrality()->GetCentralityPercentile("V0M");
888 if(fDebug>3) printf("centrality: %f\n",cent);
889 }
890 }
891
892 return GetCentralityClass(cent);
893
894}
895
896//________________________________________________________________________
897Int_t AliPWG4HighPtTrackQA::CalculateCentrality(const AliAODEvent *aod)
898{
899 //
900 // Get centrality from AOD
901 //
902
903 if(!aod) return 5;
904 Float_t cent = aod->GetHeader()->GetCentrality();
905 if(fDebug>3) printf("centrality: %f\n",cent);
906
907 return GetCentralityClass(cent);
908
909}
910
911//________________________________________________________________________
912Int_t AliPWG4HighPtTrackQA::GetCentralityClass(Float_t cent) const
913{
914 //
915 // Get centrality class
916 //
917
918 if(cent<0) return 5; // OB - cent sometimes negative
919 if(cent>80) return 4;
920 if(cent>50) return 3;
921 if(cent>30) return 2;
922 if(cent>10) return 1;
923 return 0;
924
925}
926
927//________________________________________________________________________
928void AliPWG4HighPtTrackQA::UserExec(Option_t *)
929{
930 // Main loop
931 // Called for each event
932 AliDebug(2,Form(">> AliPWG4HighPtTrackQA::UserExec \n"));
933
934 fEvent = InputEvent();
935 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
936
937 // All events without selection
938 fNEventAll->Fill(0.);
939
940 if(!SelectEvent()) {
941 // Post output data
942 PostData(1, fHistList);
943 return;
944 }
945
946
947 //Need to keep track of selected events
948 fNEventSel->Fill(0.);
949
950 fVariables = new TArrayF(fNVariables);
951
952 if(fDataType==kESD) DoAnalysisESD();
953 if(fDataType==kAOD) DoAnalysisAOD();
954
955 //Delete old fVariables
956 if(fVariables) delete fVariables;
957
958 // Post output data
959 PostData(1, fHistList);
960
961}
962
963//________________________________________________________________________
964void AliPWG4HighPtTrackQA::DoAnalysisESD()
965{
966 //
967 // Run analysis on ESD
968 //
969
970 if(!fESD) {
971 PostData(1, fHistList);
972 return;
973 }
974
975 // ---- Get MC Header information (for MC productions in pThard bins) ----
976 Double_t ptHard = 0.;
977 Double_t nTrials = 1; // trials for MC trigger weight for real data
978
979 AliMCEventHandler *eventHandlerMC = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
980 if (eventHandlerMC) {
981
982 if(eventHandlerMC->MCEvent()){
983 AliGenPythiaEventHeader* pythiaGenHeader = GetPythiaEventHeader(eventHandlerMC->MCEvent());
984 if(pythiaGenHeader){
985 nTrials = pythiaGenHeader->Trials();
986 ptHard = pythiaGenHeader->GetPtHard();
987
988 fh1PtHard->Fill(ptHard);
989 fh1PtHardTrials->Fill(ptHard,nTrials);
990
991 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
992 }
993 }
994 }
995
996 Int_t nTracks = fESD->GetNumberOfTracks();
997 AliDebug(2,Form("nTracks ESD%d", nTracks));
998
999 /*
1000 Variables to be put in fVariables
1001 0: pt
1002 1: phi
1003 2: eta
1004 3: dca2D
1005 4: dcaZ
1006 5: nClustersTPC
1007 6: nPointITS
1008 7: chi2C
1009 8: nSigmaToVertex
1010 9: trackLengthTPC
1011 10: chi2PerClusterTPC
1012 11: #crossed rows
1013 12: (#crossed rows)/(#findable clusters)
1014 13: SigmaY2
1015 14: SigmaZ2
1016 15: SigmaSnp2
1017 16: SigmaTgl2
1018 17: Sigma1Pt2
1019 18: NClustersTPCIter1
1020 19: Chi2TPCIter1
1021 20: nClustersTPCShared
1022 21: Golden Chi2 - global vs TPC constrained
1023 22: Chi2 between global and global constrained
1024 23: #crossed rows from fit map
1025 24: (#crossed rows)/(#findable clusters) from fit map
1026 25: chi2ITS
1027 */
1028
1029 for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
1030 fh1NTracksAll->Fill(0.);
1031
1032 //Get track for analysis
1033 AliESDtrack *track = 0x0;
1034 AliESDtrack *esdtrack = fESD->GetTrack(iTrack);
1035 if(!esdtrack) {
1036 fh1NTracksReject->Fill("noESDtrack",1);
1037 continue;
1038 }
1039 AliESDtrack *origtrack = new AliESDtrack(*esdtrack);
1040 if(!origtrack)
1041 continue;
1042
1043 if(fTrackType==4) {
1044 if (!(fTrackCuts->AcceptTrack(esdtrack))) {
1045 fh1NTracksReject->Fill("trackCuts",1);
1046 if(origtrack) delete origtrack;
1047 continue;
1048 }
1049 }
1050
1051 if(fTrackType==1)
1052 track = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID());
1053 else if(fTrackType==2 || fTrackType==4) {
1054 track = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(fESD),esdtrack->GetID());
1055 if(!track) {
1056 fh1NTracksReject->Fill("noTPConly",1);
1057 if(origtrack) delete origtrack;
1058 continue;
1059 }
1060 AliExternalTrackParam exParam;
1061 Bool_t relate = track->RelateToVertexTPC(fVtx,fESD->GetMagneticField(),kVeryBig,&exParam);
1062 if( !relate ) {
1063 fh1NTracksReject->Fill("relate",1);
1064 if(track) delete track;
1065 if(origtrack) delete origtrack;
1066 continue;
1067 }
1068 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
1069 }
1070 else if(fTrackType==5 || fTrackType==6) {
1071 if(fTrackCuts->AcceptTrack(esdtrack)) {
1072 if(origtrack) delete origtrack;
1073 continue;
1074 }
1075 else {
1076 if( !(fTrackCutsITSLoose->AcceptTrack(esdtrack)) && fTrackCutsTPConly->AcceptTrack(esdtrack) ) {
1077
1078 if(fTrackType==5) {
1079 //use TPConly constrained track
1080 track = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID());
1081 if(!track) {
1082 fh1NTracksReject->Fill("noTPConly",1);
1083 if(origtrack) delete origtrack;
1084 continue;
1085 }
1086 AliExternalTrackParam exParam;
1087 Bool_t relate = track->RelateToVertexTPC(fVtx,fESD->GetMagneticField(),kVeryBig,&exParam);
1088 if( !relate ) {
1089 fh1NTracksReject->Fill("relate",1);
1090 if(track) delete track;
1091 if(origtrack) delete origtrack;
1092 continue;
1093 }
1094 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
1095 }
1096 else if(fTrackType==6) {
1097 //use global constrained track
1098 track = new AliESDtrack(*esdtrack);
1099 track->Set(esdtrack->GetConstrainedParam()->GetX(),esdtrack->GetConstrainedParam()->GetAlpha(),esdtrack->GetConstrainedParam()->GetParameter(),esdtrack->GetConstrainedParam()->GetCovariance());
1100
1101 }
1102 }
1103 }
1104 }
1105 else if(fTrackType==7) {
1106 //use global constrained track
1107 track = new AliESDtrack(*esdtrack);
1108 }
1109 else
1110 track = esdtrack;
1111
1112 if(!track) {
1113 if(origtrack) delete origtrack;
1114 continue;
1115 }
1116
1117 if(fTrackType==2 || fTrackType==4 || fTrackType==5) {
1118 //Cut on chi2 of constrained fit
1119 if(track->GetConstrainedChi2TPC() > fSigmaConstrainedMax*fSigmaConstrainedMax && fSigmaConstrainedMax>0.) {
1120 fh1NTracksReject->Fill("chi2",1);
1121 if(track) delete track;
1122 if(origtrack) delete origtrack;
1123 continue;
1124 }
1125 }
1126
1127 fPtAll->Fill(track->Pt());
1128
1129 if (!(fTrackCuts->AcceptTrack(track)) && fTrackType!=4 && fTrackType!=5 && fTrackType!=6) {
1130 fh1NTracksReject->Fill("trackCuts",1);
1131 if(fTrackType==1 || fTrackType==2 || fTrackType==7) {
1132 if(track) delete track;
1133 }
1134 if(origtrack) delete origtrack;
1135 continue;
1136 }
1137
1138 if(fTrackType==7) {
1139 if(fTrackCutsITSLoose ) {
1140 if(fTrackCutsITSLoose->AcceptTrack(track) ) {
1141 if(track) delete track;
1142 if(origtrack) delete origtrack;
1143 continue;
1144 }
1145 }
1146
1147 if(esdtrack->GetConstrainedParam())
1148 track->Set(esdtrack->GetConstrainedParam()->GetX(),esdtrack->GetConstrainedParam()->GetAlpha(),esdtrack->GetConstrainedParam()->GetParameter(),esdtrack->GetConstrainedParam()->GetCovariance());
1149 }
1150
1151 if(!track) {
1152 if(fTrackType==1 || fTrackType==2 || fTrackType==4 || fTrackType==5 || fTrackType==6 || fTrackType==7) {
1153 if(track) delete track;
1154 }
1155 if(origtrack) delete origtrack;
1156 continue;
1157 }
1158
1159 fh1NTracksSel->Fill(0.);
1160
1161 fVariables->Reset(0.);
1162
1163 fVariables->SetAt(track->Pt(),0);
1164 fVariables->SetAt(track->Phi(),1);
1165 fVariables->SetAt(track->Eta(),2);
1166
1167 Float_t dca2D = 0.;
1168 Float_t dcaz = 0.;
1169
1170 if(fTrackType==1 || fTrackType==2 || fTrackType==4) {
1171 track->GetImpactParametersTPC(dca2D,dcaz); //TPConly
1172 }
1173 else
1174 track->GetImpactParameters(dca2D,dcaz); //Global
1175
1176 fVariables->SetAt(dca2D,3);
1177 fVariables->SetAt(dcaz,4);
1178
1179 fVariables->SetAt((float)track->GetTPCNcls(),5);
1180
1181 Int_t nPointITS = 0;
1182 fITSClusterMap = track->GetITSClusterMap();
1183 UChar_t itsMap = track->GetITSClusterMap();
1184 for (Int_t i=0; i < 6; i++) {
1185 if (itsMap & (1 << i))
1186 nPointITS ++;
1187 }
1188 fVariables->SetAt((float)nPointITS,6);
1189 Float_t chi2C = (float)track->GetConstrainedChi2();
1190 if(fTrackType==1 || fTrackType==2 || fTrackType==4)
1191 chi2C = (float)track->GetConstrainedChi2TPC();
1192 fVariables->SetAt(chi2C,7);
1193 fVariables->SetAt(fTrackCuts->GetSigmaToVertex(track),8);// Calculates the number of sigma to the vertex for a track.
1194
1195 fVariables->SetAt(GetTrackLengthTPC(track),9);
1196
1197 if(fVariables->At(5)>0.) fVariables->SetAt(track->GetTPCchi2()/fVariables->At(5),10);
1198
1199 //fVariables->SetAt(track->GetTPCClusterInfo(2,1),11); //#crossed rows
1200 fVariables->SetAt(track->GetTPCCrossedRows(),11); //#crossed rows
1201
1202 Float_t crossedRowsTPCNClsF = 1.;//track->GetTPCClusterInfo(2,0);
1203 if(track->GetTPCNclsF()>0.) crossedRowsTPCNClsF = fVariables->At(11)/track->GetTPCNclsF();
1204 fVariables->SetAt(crossedRowsTPCNClsF,12);//(#crossed rows)/(#findable clusters)
1205 fVariables->SetAt(track->GetSigmaY2(),13);
1206 fVariables->SetAt(track->GetSigmaZ2(),14);
1207 fVariables->SetAt(track->GetSigmaSnp2(),15);
1208 fVariables->SetAt(track->GetSigmaTgl2(),16);
1209 fVariables->SetAt(track->GetSigma1Pt2(),17);
1210
1211 fVariables->SetAt(track->GetTPCNclsIter1(),18);
1212 fVariables->SetAt(track->GetTPCchi2Iter1(),19);
1213
1214 fVariables->SetAt(track->GetTPCnclsS(),20);
1215
1216 Float_t chi2Gold = origtrack->GetChi2TPCConstrainedVsGlobal(fVtx);//GetGoldenChi2(origtrack);
1217 Float_t chi2GGC = GetGGCChi2(origtrack);
1218
1219 fVariables->SetAt(chi2Gold,21);
1220 fVariables->SetAt(chi2GGC,22);
1221
1222 fVariables->SetAt(GetTPCClusterInfoFitMap(track,2,1),23);
1223 Float_t crossedRowsTPCNClsFFit = 1.;
1224 if(track->GetTPCNclsF()>0.) crossedRowsTPCNClsFFit = fVariables->At(23)/track->GetTPCNclsF();
1225 fVariables->SetAt(crossedRowsTPCNClsFFit,24);
1226
1227 fVariables->SetAt(track->GetITSchi2(),25);
1228
1229 TBits fitmap = track->GetTPCFitMap();
1230 fPtNClustersNClustersFitMap->Fill(track->Pt(),track->GetTPCNcls(),(float)fitmap.CountBits());
1231
1232 FillHistograms();
1233
1234 // int mult = fTrackCuts->CountAcceptedTracks(fESD);
1235
1236 if(fTrackType==1 || fTrackType==2 || fTrackType==4 || fTrackType==5 || fTrackType==6 || fTrackType==7) {
1237 if(track) delete track;
1238 }
1239 if(origtrack) delete origtrack;
1240
1241 }//track loop
1242
1243}
1244
1245//________________________________________________________________________
1246void AliPWG4HighPtTrackQA::DoAnalysisAOD()
1247{
1248 //
1249 // Do QA on AOD input
1250 //
1251 AliAODEvent *aod = dynamic_cast<AliAODEvent*>(fEvent);
1252 if(!aod) return;
1253 AliExternalTrackParam exParam;
1254 for (Int_t iTrack = 0; iTrack < fEvent->GetNumberOfTracks(); iTrack++) {
1255
1256 AliAODTrack *aodtrack = aod->GetTrack(iTrack);
1257 if( !aodtrack->TestFilterMask(fFilterMask) ) {
1258 fh1NTracksReject->Fill("noHybridTrack",1);
1259 continue;
1260 }
1261
1262 if(!fIncludeNoITS) {
1263 if ((aodtrack->GetStatus()&AliESDtrack::kITSrefit)==0) {
1264 fh1NTracksReject->Fill("noITSrefit",1);
1265 continue;
1266 }
1267 }
1268
1269 fVariables->Reset(0.);
1270
1271 fVariables->SetAt(aodtrack->Pt(),0);
1272 fVariables->SetAt(aodtrack->Phi(),1);
1273 fVariables->SetAt(aodtrack->Eta(),2);
1274
1275 Double_t dca[2] = {0.,0.};
1276 if(aodtrack->IsGlobalConstrained()) {
1277 dca[0] = aodtrack->DCA();
1278 dca[1] = aodtrack->ZAtDCA();
1279 } else {
1280 Double_t v[3] = {0};
1281 Double_t pos[3] = {0};
1282 fVtxAOD->GetXYZ(v);
1283 aodtrack->GetXYZ(pos);
1284 dca[0] = pos[0] - v[0];
1285 dca[1] = pos[1] - v[1];
1286 }
1287 fVariables->SetAt(dca[0],3);
1288 fVariables->SetAt(dca[1],4);
1289 fVariables->SetAt((float)aodtrack->GetTPCNcls(),5);
1290 fVariables->SetAt((float)aodtrack->GetITSNcls(),6);
1291 fVariables->SetAt(0.,7); //ConstrainedChi2TPC -> not available in AOD
1292 fVariables->SetAt(0.,8);
1293 fVariables->SetAt(GetTrackLengthTPC(aodtrack),9);
1294 Float_t chi2pndf = aodtrack->Chi2perNDF();
1295 //if(fVariables->At(5)>0.) chi2pndf = aodtrack->GetTPCchi2()/fVariables->At(5);
1296 fVariables->SetAt(chi2pndf,10);
1297 fVariables->SetAt(GetTPCClusterInfo(aodtrack,2,1,0,159,kFALSE),11);
1298 Float_t crossedRowsTPCNClsF = 0.;
1299 if(aodtrack->GetTPCNclsF()>0.) crossedRowsTPCNClsF = fVariables->At(11)/aodtrack->GetTPCNclsF();
1300 fVariables->SetAt(crossedRowsTPCNClsF,12);
1301
1302 //get covariance matrix
1303 Double_t cov[21] = {0,};
1304 aodtrack->GetCovMatrix(cov);
1305 Double_t pxpypz[3] = {0,};
1306 aodtrack->PxPyPz(pxpypz);
1307 Double_t xyz[3] = {0,};
1308 aodtrack->GetXYZ(xyz);
1309 Short_t sign = aodtrack->Charge();
1310 exParam.Set(xyz,pxpypz,cov,sign);
1311
1312 fVariables->SetAt(exParam.GetSigmaY2(),13);
1313 fVariables->SetAt(exParam.GetSigmaZ2(),14);
1314 fVariables->SetAt(exParam.GetSigmaSnp2(),15);
1315 fVariables->SetAt(exParam.GetSigmaTgl2(),16);
1316 fVariables->SetAt(exParam.GetSigma1Pt2(),17);
1317
1318 fVariables->SetAt(0.,18); //NClustersTPCIter1
1319 fVariables->SetAt(0.,19); //Chi2TPCIter1
1320
1321 TBits sharedClusterMap = aodtrack->GetTPCSharedMap();
1322 fVariables->SetAt(sharedClusterMap.CountBits(),20);
1323
1324 fVariables->SetAt(0.,21); //not available in AOD golden chi2
1325 fVariables->SetAt(0.,22); //not available in AOD Chi2 between global and global constrained
1326
1327 fVariables->SetAt(GetTPCClusterInfo(aodtrack,2,1,0,159,kTRUE),23); //not available in AOD #crossed rows from fit map
1328 Float_t crossedRowsTPCNClsFFit = 0.;
1329 if(aodtrack->GetTPCNclsF()>0.) crossedRowsTPCNClsFFit = fVariables->At(23)/aodtrack->GetTPCNclsF();
1330 fVariables->SetAt(crossedRowsTPCNClsFFit,24); //(#crossed rows)/(#findable clusters) from fit map
1331
1332 fVariables->SetAt(0.,25);
1333
1334 fPtAll->Fill(fVariables->At(0));
1335
1336 FillHistograms();
1337 }
1338}
1339
1340//________________________________________________________________________
1341void AliPWG4HighPtTrackQA::FillHistograms()
1342{
1343 //
1344 // Fill all QA histograms
1345 //
1346
1347 fPtSel->Fill(fVariables->At(0));
1348 fPtPhi->Fill(fVariables->At(0),fVariables->At(1));
1349 fPtEta->Fill(fVariables->At(0),fVariables->At(2));
1350 fPtEtaPhi->Fill(fVariables->At(0),fVariables->At(2),fVariables->At(1));
1351 fPtDCA2D->Fill(fVariables->At(0),fVariables->At(3));
1352 fPtDCAZ->Fill(fVariables->At(0),fVariables->At(4));
1353 fPtNClustersTPC->Fill(fVariables->At(0),fVariables->At(5));
1354 fPtNClustersTPCPhi->Fill(fVariables->At(1),fVariables->At(5));
1355 fPtNPointITS->Fill(fVariables->At(0),fVariables->At(6));
1356 fPtNPointITSPhi->Fill(fVariables->At(0),fVariables->At(6),fVariables->At(1));
1357
1358 fPtNClustersTPCIter1->Fill(fVariables->At(0),fVariables->At(18));
1359 fPtNClustersTPCIter1Phi->Fill(fVariables->At(0),fVariables->At(18),fVariables->At(1));
1360 fPtNClustersTPCShared->Fill(fVariables->At(0),fVariables->At(20));
1361 if(fVariables->At(5)>0.)
1362 fPtNClustersTPCSharedFrac->Fill(fVariables->At(0),fVariables->At(20)/fVariables->At(5));
1363
1364 if(fVariables->At(18)>0.)
1365 fPtChi2PerClusterTPCIter1->Fill(fVariables->At(0),fVariables->At(19)/fVariables->At(18));
1366
1367 fPtChi2C->Fill(fVariables->At(0),fVariables->At(7));
1368 fPtNSigmaToVertex->Fill(fVariables->At(0),fVariables->At(8));
1369 fPtRelUncertainty1Pt->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)));
1370 fPtRelUncertainty1PtNClus->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)),fVariables->At(5));
1371 fPtRelUncertainty1PtNClusIter1->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)),fVariables->At(18));
1372 fPtRelUncertainty1PtNPointITS->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)),fVariables->At(6));
1373
1374 fPtRelUncertainty1PtITSClusterMap->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)),(int)fITSClusterMap);
1375
1376 fPtRelUncertainty1PtChi2->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)),fVariables->At(10));
1377 if(fVariables->At(18)>0.)
1378 fPtRelUncertainty1PtChi2Iter1->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)),fVariables->At(19)/fVariables->At(18));
1379 fPtRelUncertainty1PtPhi->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)),fVariables->At(1));
1380
1381 fPtSigmaY2->Fill(1./fVariables->At(0),TMath::Sqrt(fVariables->At(13)));
1382 fPtSigmaZ2->Fill(1./fVariables->At(0),TMath::Sqrt(fVariables->At(14)));
1383 fPtSigmaSnp2->Fill(1./fVariables->At(0),TMath::Sqrt(fVariables->At(15)));
1384 fPtSigmaTgl2->Fill(1./fVariables->At(0),TMath::Sqrt(fVariables->At(16)));
1385 fPtSigma1Pt2->Fill(1./fVariables->At(0),TMath::Sqrt(fVariables->At(17)));
1386
1387 fProfPtSigmaY2->Fill(1./fVariables->At(0),TMath::Sqrt(fVariables->At(13)));
1388 fProfPtSigmaZ2->Fill(1./fVariables->At(0),TMath::Sqrt(fVariables->At(14)));
1389 fProfPtSigmaSnp2->Fill(1./fVariables->At(0),TMath::Sqrt(fVariables->At(15)));
1390 fProfPtSigmaTgl2->Fill(1./fVariables->At(0),TMath::Sqrt(fVariables->At(16)));
1391 fProfPtSigma1Pt2->Fill(1./fVariables->At(0),TMath::Sqrt(fVariables->At(17)));
1392 fProfPtSigma1Pt->Fill(fVariables->At(0),TMath::Sqrt(fVariables->At(17)));
1393 fProfPtPtSigma1Pt->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)));
1394
1395 fPtChi2PerClusterTPC->Fill(fVariables->At(0),fVariables->At(10));
1396 fPtNCrossedRows->Fill(fVariables->At(0),fVariables->At(11));
1397 fPtNCrossedRowsPhi->Fill(fVariables->At(0),fVariables->At(11),fVariables->At(1));
1398 fPtNCrossedRowsNClusFPhi->Fill(fVariables->At(0),fVariables->At(12),fVariables->At(1));
1399 fPtNCrRNCrRNClusF->Fill(fVariables->At(0),fVariables->At(11),fVariables->At(12));
1400
1401 fPtChi2Gold->Fill(fVariables->At(0),fVariables->At(21));
1402 fPtChi2GGC->Fill(fVariables->At(0),fVariables->At(22));
1403
1404 fPtChi2GoldPhi->Fill(fVariables->At(0),fVariables->At(21),fVariables->At(1));
1405 fPtChi2GGCPhi->Fill(fVariables->At(0),fVariables->At(22),fVariables->At(1));
1406
1407 fChi2GoldChi2GGC->Fill(fVariables->At(21),fVariables->At(22));
1408
1409 fPtNCrossedRowsFit->Fill(fVariables->At(0),fVariables->At(23));
1410 fPtNCrossedRowsFitPhi->Fill(fVariables->At(0),fVariables->At(23),fVariables->At(1));
1411 fPtNCrossedRowsNClusFFitPhi->Fill(fVariables->At(0),fVariables->At(24),fVariables->At(1));
1412 fNCrossedRowsNCrossedRowsFit->Fill(fVariables->At(11),fVariables->At(23));
1413
1414 fNClustersNCrossedRows->Fill(fVariables->At(5),fVariables->At(11));
1415 fNClustersNCrossedRowsFit->Fill(fVariables->At(5),fVariables->At(23));
1416
1417 fPtRelUncertainty1PtNCrossedRows->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)),fVariables->At(11));
1418 fPtRelUncertainty1PtNCrossedRowsFit->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)),fVariables->At(23));
1419
1420 if(fVariables->At(6)>0.)
1421 fPtChi2ITSPhi->Fill(fVariables->At(0),fVariables->At(25)/fVariables->At(6),fVariables->At(1));
1422
1423}
1424
1425//________________________________________________________________________
1426Bool_t AliPWG4HighPtTrackQA::PythiaInfoFromFile(const char* currFile,Float_t &fXsec,Float_t &fTrials)
1427{
1428 //
1429 // get the cross section and the trails either from pyxsec.root or from pysec_hists.root
1430 // This is to called in Notify and should provide the path to the AOD/ESD file
1431 // Copied from AliAnalysisTaskJetSpectrum2
1432 //
1433
1434 TString file(currFile);
1435 fXsec = 0;
1436 fTrials = 1;
1437
1438 if(file.Contains("root_archive.zip#")){
1439 Ssiz_t pos1 = file.Index("root_archive",12,TString::kExact);
1440 Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
1441 file.Replace(pos+1,20,"");
1442 }
1443 else {
1444 // not an archive take the basename....
1445 file.ReplaceAll(gSystem->BaseName(file.Data()),"");
1446 }
1447
1448 TFile *fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec.root")); // problem that we cannot really test the existance of a file in a archive so we have to lvie with open error message from root
1449 if(!fxsec){
1450 // next trial fetch the histgram file
1451 fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
1452 if(!fxsec){
1453 // not a severe condition but inciate that we have no information
1454 return kFALSE;
1455 }
1456 else{
1457 // find the tlist we want to be independtent of the name so use the Tkey
1458 TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
1459 if(!key){
1460 fxsec->Close();
1461 return kFALSE;
1462 }
1463 TList *list = dynamic_cast<TList*>(key->ReadObj());
1464 if(!list){
1465 fxsec->Close();
1466 return kFALSE;
1467 }
1468 fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
1469 fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
1470 fxsec->Close();
1471 }
1472 } // no tree pyxsec.root
1473 else {
1474 TTree *xtree = (TTree*)fxsec->Get("Xsection");
1475 if(!xtree){
1476 fxsec->Close();
1477 return kFALSE;
1478 }
1479 UInt_t ntrials = 0;
1480 Double_t xsection = 0;
1481 xtree->SetBranchAddress("xsection",&xsection);
1482 xtree->SetBranchAddress("ntrials",&ntrials);
1483 xtree->GetEntry(0);
1484 fTrials = ntrials;
1485 fXsec = xsection;
1486 fxsec->Close();
1487 }
1488 return kTRUE;
1489}
1490
1491//________________________________________________________________________
1492Bool_t AliPWG4HighPtTrackQA::Notify()
1493{
1494 //
1495 // Implemented Notify() to read the cross sections
1496 // and number of trials from pyxsec.root
1497 // Copied from AliAnalysisTaskJetSpectrum2
1498 //
1499
1500 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
1501 Float_t xsection = 0;
1502 Float_t ftrials = 1;
1503
1504 fAvgTrials = 1;
1505 if(tree){
1506 TFile *curfile = tree->GetCurrentFile();
1507 if (!curfile) {
1508 Error("Notify","No current file");
1509 return kFALSE;
1510 }
1511 if(!fh1Xsec||!fh1Trials){
1512 // Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
1513 return kFALSE;
1514 }
1515 PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
1516 fh1Xsec->Fill("<#sigma>",xsection);
1517 // construct a poor man average trials
1518 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
1519 if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
1520 }
1521 return kTRUE;
1522}
1523
1524//________________________________________________________________________
1525AliGenPythiaEventHeader* AliPWG4HighPtTrackQA::GetPythiaEventHeader(const AliMCEvent *mcEvent)
1526{
1527
1528 if(!mcEvent)return 0;
1529 AliGenEventHeader* genHeader = mcEvent->GenEventHeader();
1530 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
1531 if(!pythiaGenHeader){
1532 // cocktail ??
1533 AliGenCocktailEventHeader* genCocktailHeader = dynamic_cast<AliGenCocktailEventHeader*>(genHeader);
1534
1535 if (!genCocktailHeader) {
1536 // AliWarningGeneral(Form(" %s:%d",(char*)__FILE__,__LINE__),"Unknown header type (not Pythia or Cocktail)");
1537 // AliWarning(Form("%s %d: Unknown header type (not Pythia or Cocktail)",(char*)__FILE__,__LINE__));
1538 return 0;
1539 }
1540 TList* headerList = genCocktailHeader->GetHeaders();
1541 for (Int_t i=0; i<headerList->GetEntries(); i++) {
1542 pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(headerList->At(i));
1543 if (pythiaGenHeader)
1544 break;
1545 }
1546 if(!pythiaGenHeader){
1547 AliWarningGeneral(Form(" %s:%d",(char*)__FILE__,__LINE__),"Pythia event header not found");
1548 return 0;
1549 }
1550 }
1551 return pythiaGenHeader;
1552
1553}
1554
1555//_______________________________________________________________________
1556Float_t AliPWG4HighPtTrackQA::GetTPCClusterInfo(const AliAODTrack *tr,Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1, Bool_t useFitMap) const
1557{
1558 //MV: copied from AliESDtrack since method is not available in AliAODTrack
1559
1560 //
1561 // TPC cluster information
1562 // type 0: get fraction of found/findable clusters with neighbourhood definition
1563 // 1: findable clusters with neighbourhood definition
1564 // 2: found clusters
1565 //
1566 // definition of findable clusters:
1567 // a cluster is defined as findable if there is another cluster
1568 // within +- nNeighbours pad rows. The idea is to overcome threshold
1569 // effects with a very simple algorithm.
1570 //
1571
1572 TBits fTPCClusterMap = 0;
1573 if(useFitMap)
1574 fTPCClusterMap = tr->GetTPCFitMap();
1575 else
1576 fTPCClusterMap = tr->GetTPCClusterMap();
1577
1578 if (type==2) return fTPCClusterMap.CountBits();
1579
1580 Int_t found=0;
1581 Int_t findable=0;
1582 Int_t last=-nNeighbours;
1583
1584 for (Int_t i=row0; i<row1; ++i){
1585 //look to current row
1586 if (fTPCClusterMap[i]) {
1587 last=i;
1588 ++found;
1589 ++findable;
1590 continue;
1591 }
1592 //look to nNeighbours before
1593 if ((i-last)<=nNeighbours) {
1594 ++findable;
1595 continue;
1596 }
1597 //look to nNeighbours after
1598 for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
1599 if (fTPCClusterMap[j]){
1600 ++findable;
1601 break;
1602 }
1603 }
1604 }
1605 if (type==1) return findable;
1606
1607 if (type==0){
1608 Float_t fraction=0;
1609 if (findable>0)
1610 fraction=(Float_t)found/(Float_t)findable;
1611 else
1612 fraction=0;
1613 return fraction;
1614 }
1615 return 0; // undefined type - default value
1616}
1617
1618//_______________________________________________________________________
1619Float_t AliPWG4HighPtTrackQA::GetTPCClusterInfoFitMap(const AliESDtrack *tr,Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1) const
1620{
1621 //
1622 // TPC cluster information from fit map
1623 // type 0: get fraction of found/findable clusters with neighbourhood definition
1624 // 1: findable clusters with neighbourhood definition
1625 // 2: found clusters
1626 //
1627 // definition of findable clusters:
1628 // a cluster is defined as findable if there is another cluster
1629 // within +- nNeighbours pad rows. The idea is to overcome threshold
1630 // effects with a very simple algorithm.
1631 //
1632
1633 TBits fTPCFitMap = tr->GetTPCFitMap();
1634 if (type==2) return fTPCFitMap.CountBits();
1635
1636 Int_t found=0;
1637 Int_t findable=0;
1638 Int_t last=-nNeighbours;
1639
1640 for (Int_t i=row0; i<row1; ++i){
1641 //look to current row
1642 if (fTPCFitMap[i]) {
1643 last=i;
1644 ++found;
1645 ++findable;
1646 continue;
1647 }
1648 //look to nNeighbours before
1649 if ((i-last)<=nNeighbours) {
1650 ++findable;
1651 continue;
1652 }
1653 //look to nNeighbours after
1654 for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
1655 if (fTPCFitMap[j]){
1656 ++findable;
1657 break;
1658 }
1659 }
1660 }
1661 if (type==1) return findable;
1662
1663 if (type==0){
1664 Float_t fraction=0;
1665 if (findable>0)
1666 fraction=(Float_t)found/(Float_t)findable;
1667 else
1668 fraction=0;
1669 return fraction;
1670 }
1671 return 0; // undefined type - default value
1672}
1673
1674//_______________________________________________________________________
1675Int_t AliPWG4HighPtTrackQA::GetTrackLengthTPC(const AliESDtrack *track) const
1676{
1677 //
1678 // returns distance between 1st and last hit in TPC
1679 // distance given in number of padrows
1680 //
1681
1682 TBits fTPCClusterMap = track->GetTPCClusterMap();
1683 int firstHit = 0;
1684 int lastHit = 0;
1685
1686 for(int i=0; i<=159; i++) {
1687 if(fTPCClusterMap[i]>0) firstHit = i;
1688 }
1689 for(int i=159; i>=0; i--) {
1690 if(fTPCClusterMap[i]>0) lastHit = i;
1691 }
1692
1693 Int_t trackLength = lastHit - firstHit;
1694
1695 return trackLength;
1696}
1697
1698//_______________________________________________________________________
1699Int_t AliPWG4HighPtTrackQA::GetTrackLengthTPC(const AliAODTrack *track) const
1700{
1701 //
1702 // returns distance between 1st and last hit in TPC
1703 // distance given in number of padrows
1704 //
1705
1706 TBits fTPCClusterMap = track->GetTPCClusterMap();
1707 int firstHit = 0;
1708 int lastHit = 0;
1709
1710 for(int i=0; i<=159; i++) {
1711 if(fTPCClusterMap[i]>0) firstHit = i;
1712 }
1713 for(int i=159; i>=0; i--) {
1714 if(fTPCClusterMap[i]>0) lastHit = i;
1715 }
1716
1717 Int_t trackLength = lastHit - firstHit;
1718
1719 return trackLength;
1720}
1721
1722//_______________________________________________________________________
1723Float_t AliPWG4HighPtTrackQA::GetGoldenChi2(AliESDtrack *origtrack)
1724{
1725 //
1726 // Return chi2 between global and TPC constrained track
1727 // track should be the global unconstrained track
1728 //
1729
1730 Float_t chi2Gold = 0.;
1731
1732 AliESDtrack *tpcTrack = 0x0;
1733 tpcTrack = AliESDtrackCuts::GetTPCOnlyTrack(fESD,origtrack->GetID());
1734 if(tpcTrack) {
1735 AliExternalTrackParam exParam;
1736 Bool_t relate = tpcTrack->RelateToVertexTPC(fVtx,fESD->GetMagneticField(),kVeryBig,&exParam);
1737 if( relate ) {
1738 tpcTrack->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
1739 // Double_t pTPC[2],covTPC[3]; tpcTrack->PropagateToDCA(fVtx, fESD->GetMagneticField(), 10000, pTPC, covTPC);
1740 }
1741
1742 tpcTrack->Propagate(origtrack->GetAlpha(), origtrack->GetX(), fESD->GetMagneticField());
1743 chi2Gold = (Float_t)origtrack->GetPredictedChi2(tpcTrack);
1744 }
1745
1746 if(tpcTrack) delete tpcTrack;
1747
1748 return chi2Gold;
1749
1750}
1751
1752//_______________________________________________________________________
1753Float_t AliPWG4HighPtTrackQA::GetGGCChi2(AliESDtrack *origtrack)
1754{
1755 //
1756 // Return chi2 between global and global constrained track
1757 // track should be the global unconstrained track
1758 //
1759
1760 Float_t chi2GGC = 0.;
1761
1762 AliESDtrack *esdtrackC = new AliESDtrack(*origtrack);
1763 if(esdtrackC) {
1764 if(origtrack->GetConstrainedParam()) {
1765 esdtrackC->Set(origtrack->GetConstrainedParam()->GetX(),origtrack->GetConstrainedParam()->GetAlpha(),origtrack->GetConstrainedParam()->GetParameter(),origtrack->GetConstrainedParam()->GetCovariance());
1766 chi2GGC = (Float_t)origtrack->GetPredictedChi2(esdtrackC);
1767 }
1768 delete esdtrackC;
1769 }
1770
1771 return chi2GGC;
1772
1773}
1774
1775//________________________________________________________________________
1776void AliPWG4HighPtTrackQA::Terminate(Option_t *)
1777{
1778 // The Terminate() function is the last function to be called during
1779 // a query. It always runs on the client, it can be used to present
1780 // the results graphically or save the results to file.
1781
1782}
1783
1784#endif