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