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