1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
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 //-----------------------------------------------------------------------
24 #ifndef ALIPWG4HIGHPTTRACKQA_CXX
25 #define ALIPWG4HIGHPTTRACKQA_CXX
27 #include "AliPWG4HighPtTrackQA.h"
41 #include "AliAnalysisManager.h"
42 #include "AliESDInputHandler.h"
43 #include "AliMCEvent.h"
44 #include "AliMCEventHandler.h"
46 #include "AliESDtrack.h"
47 #include "AliESDtrackCuts.h"
48 #include "AliExternalTrackParam.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"
57 using namespace std; //required for resolving the 'cout' symbol
59 ClassImp(AliPWG4HighPtTrackQA)
61 AliPWG4HighPtTrackQA::AliPWG4HighPtTrackQA()
62 : AliAnalysisTaskSE(),
68 fTrackCutsITSLoose(0x0),
69 fTrackCutsTPConly(0x0),
72 fSigmaConstrainedMax(-1.),
89 fh1NTracksReject(0x0),
98 fPtNClustersTPCIter1(0x0),
99 fPtNClustersTPCIter1Phi(0x0),
100 fPtNClustersTPCShared(0x0),
101 fPtNClustersTPCSharedFrac(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),
132 fChi2GoldChi2GGC(0x0),
140 fProfPtSigmaSnp2(0x0),
141 fProfPtSigmaTgl2(0x0),
142 fProfPtSigma1Pt2(0x0),
143 fProfPtSigma1Pt(0x0),
144 fProfPtPtSigma1Pt(0x0),
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.;
160 //________________________________________________________________________
161 AliPWG4HighPtTrackQA::AliPWG4HighPtTrackQA(const char *name):
162 AliAnalysisTaskSE(name),
168 fTrackCutsITSLoose(0x0),
169 fTrackCutsTPConly(0x0),
172 fSigmaConstrainedMax(-1.),
189 fh1NTracksReject(0x0),
197 fPtNClustersTPC(0x0),
198 fPtNClustersTPCIter1(0x0),
199 fPtNClustersTPCIter1Phi(0x0),
200 fPtNClustersTPCShared(0x0),
201 fPtNClustersTPCSharedFrac(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),
232 fChi2GoldChi2GGC(0x0),
240 fProfPtSigmaSnp2(0x0),
241 fProfPtSigmaTgl2(0x0),
242 fProfPtSigma1Pt2(0x0),
243 fProfPtSigma1Pt(0x0),
244 fProfPtPtSigma1Pt(0x0),
248 // Constructor. Initialization of Inputs and Outputs
250 AliDebug(2,Form("AliPWG4HighPtTrackQA Calling Constructor"));
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.;
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());
267 //________________________________________________________________________
268 void AliPWG4HighPtTrackQA::SetPtBinEdges(Int_t region, Double_t ptmax, Double_t ptBinWidth) {
270 // Set variable bin sizes for pT axis in histos
274 fPtBinEdges[region][0] = ptmax;
275 fPtBinEdges[region][1] = ptBinWidth;
278 AliError("Only 3 regions alowed. Use region 0/1/2\n");
284 //________________________________________________________________________
285 void AliPWG4HighPtTrackQA::UserCreateOutputObjects() {
286 //Create output objects
287 AliDebug(2,Form(">> AliPWG4HighPtTrackQA::UserCreateOutputObjects \n"));
289 Bool_t oldStatus = TH1::AddDirectoryStatus();
290 TH1::AddDirectory(kFALSE);
293 fHistList = new TList();
294 fHistList->SetOwner(kTRUE);
296 Float_t fgkPtMin = 0.;
297 // Float_t fgkPtMax = fPtMax;
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) ;
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 ;
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 ;
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 ;
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) {
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 ;
346 Int_t fgkNDCAZBins=80;
347 Float_t fgkDCAZMin = -2.;
348 Float_t fgkDCAZMax = 2.;
349 if(fTrackType==1 || fTrackType==2 || fTrackType==4) {
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 ;
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 ;
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 ;
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 ;
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;
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;
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);
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 ;
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 ;
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 ;
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];
428 for(Int_t i=0; i<=fgkN1PtBins; i++) {
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);
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 ;
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 ;
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 ;
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 ;
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 ;
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);
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);
489 fh1Centrality = new TH1F("fh1Centrality","fh1Centrality; Centrality %",100,0,100);
490 fHistList->Add(fh1Centrality);
492 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
493 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
494 fHistList->Add(fh1Xsec);
496 fh1Trials = new TH1F("fh1Trials","trials root file",1,0,1);
497 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
498 fHistList->Add(fh1Trials);
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);
505 fh1NTracksAll = new TH1F("fh1NTracksAll","fh1NTracksAll",1,-0.5,0.5);
506 fHistList->Add(fh1NTracksAll);
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);
520 fh1NTracksSel = new TH1F("fh1NTracksSel","fh1NTracksSel",1,-0.5,0.5);
521 fHistList->Add(fh1NTracksSel);
523 fPtAll = new TH1F("fPtAll","PtAll",fgkNPtBins, binsPt);
524 fHistList->Add(fPtAll);
525 fPtSel = new TH1F("fPtSel","PtSel",fgkNPtBins, binsPt);
526 fHistList->Add(fPtSel);
528 fPtPhi = new TH2F("fPtPhi","fPtPhi",fgkNPtBins,binsPt,fgkNPhiBins,binsPhi);
529 fHistList->Add(fPtPhi);
531 fPtEta = new TH2F("fPtEta","fPtEta",fgkNPtBins,binsPt,fgkNEtaBins,binsEta);
532 fHistList->Add(fPtEta);
534 fPtDCA2D = new TH2F("fPtDCA2D","fPtDCA2D",fgkNPtBins,binsPt,fgkNDCA2DBins,binsDCA2D);
535 fHistList->Add(fPtDCA2D);
537 fPtDCAZ = new TH2F("fPtDCAZ","fPtDCAZ",fgkNPtBins,binsPt,fgkNDCAZBins,binsDCAZ);
538 fHistList->Add(fPtDCAZ);
540 fPtNClustersTPC = new TH2F("fPtNClustersTPC","fPtNClustersTPC",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC);
541 fHistList->Add(fPtNClustersTPC);
543 fPtNClustersTPCIter1 = new TH2F("fPtNClustersTPCIter1","fPtNClustersTPCIter1",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC);
544 fHistList->Add(fPtNClustersTPCIter1);
546 fPtNClustersTPCIter1Phi = new TH3F("fPtNClustersTPCIter1Phi","fPtNClustersTPCIter1Phi",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC,fgkNPhiBins,binsPhi);
547 fHistList->Add(fPtNClustersTPCIter1Phi);
549 fPtNClustersTPCShared = new TH2F("fPtNClustersTPCShared","fPtNClustersTPCShared",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC);
550 fHistList->Add(fPtNClustersTPCShared);
552 fPtNClustersTPCSharedFrac = new TH2F("fPtNClustersTPCSharedFrac","fPtNClustersTPCSharedFrac",fgkNPtBins,binsPt,fgkNSigma1Pt2Bins,binsSigma1Pt2);
553 fHistList->Add(fPtNClustersTPCSharedFrac);
555 fPtNPointITS = new TH2F("fPtNPointITS","fPtNPointITS",fgkNPtBins,binsPt,fgkNNPointITSBins,binsNPointITS);
556 fHistList->Add(fPtNPointITS);
558 fPtChi2C = new TH2F("fPtChi2C","fPtChi2C",fgkNPtBins,binsPt,fgkNChi2CBins,binsChi2C);
559 fHistList->Add(fPtChi2C);
561 fPtNSigmaToVertex = new TH2F("fPtNSigmaToVertex","fPtNSigmaToVertex",fgkNPtBins,binsPt,fgkNNSigmaToVertexBins,binsNSigmaToVertex);
562 fHistList->Add(fPtNSigmaToVertex);
564 fPtRelUncertainty1Pt = new TH2F("fPtRelUncertainty1Pt","fPtRelUncertainty1Pt",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty);
565 fHistList->Add(fPtRelUncertainty1Pt);
567 fPtRelUncertainty1PtNClus = new TH3F("fPtRelUncertainty1PtNClus","fPtRelUncertainty1PtNClus",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNNClustersTPCBins,binsNClustersTPC);
568 fHistList->Add(fPtRelUncertainty1PtNClus);
570 fPtRelUncertainty1PtNClusIter1 = new TH3F("fPtRelUncertainty1PtNClusIter1","fPtRelUncertainty1PtNClusIter1",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNNClustersTPCBins,binsNClustersTPC);
571 fHistList->Add(fPtRelUncertainty1PtNClusIter1);
573 fPtRelUncertainty1PtNPointITS = new TH3F("fPtRelUncertainty1PtNPointITS","fPtRelUncertainty1PtNPointITS",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNNPointITSBins,binsNPointITS);
574 fHistList->Add(fPtRelUncertainty1PtNPointITS);
576 fPtRelUncertainty1PtITSClusterMap = new TH3F("fPtRelUncertainty1PtITSClusterMap","fPtRelUncertainty1PtITSClusterMap",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNITSClusterMapBins,binsITSClusterMap);
577 fHistList->Add(fPtRelUncertainty1PtITSClusterMap);
579 fPtRelUncertainty1PtChi2 = new TH3F("fPtRelUncertainty1PtChi2","fPtRelUncertainty1PtChi2",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNChi2PerClusBins,binsChi2PerClus);
580 fHistList->Add(fPtRelUncertainty1PtChi2);
582 fPtRelUncertainty1PtChi2Iter1 = new TH3F("fPtRelUncertainty1PtChi2Iter1","fPtRelUncertainty1PtChi2Iter1",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNChi2PerClusBins,binsChi2PerClus);
583 fHistList->Add(fPtRelUncertainty1PtChi2Iter1);
585 fPtRelUncertainty1PtPhi = new TH3F("fPtRelUncertainty1PtPhi","fPtRelUncertainty1PtPhi",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNPhiBins,binsPhi);
586 fHistList->Add(fPtRelUncertainty1PtPhi);
588 fPtUncertainty1Pt = new TH2F("fPtUncertainty1Pt","fPtUncertainty1Pt",fgkNPtBins,binsPt,fgkNUncertainty1PtBins,binsUncertainty1Pt);
589 fHistList->Add(fPtUncertainty1Pt);
591 fPtChi2PerClusterTPC = new TH2F("fPtChi2PerClusterTPC","fPtChi2PerClusterTPC",fgkNPtBins,binsPt,fgkNChi2PerClusBins,binsChi2PerClus);
592 fHistList->Add(fPtChi2PerClusterTPC);
594 fPtChi2PerClusterTPCIter1 = new TH2F("fPtChi2PerClusterTPCIter1","fPtChi2PerClusterTPCIter1",fgkNPtBins,binsPt,fgkNChi2PerClusBins,binsChi2PerClus);
595 fHistList->Add(fPtChi2PerClusterTPCIter1);
597 fPtNCrossedRows = new TH2F("fPtNCrossedRows","fPtNCrossedRows",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC);
598 fHistList->Add(fPtNCrossedRows);
600 fPtNCrossedRowsPhi = new TH3F("fPtNCrossedRowsPhi","fPtNCrossedRowsPhi",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC,fgkNPhiBins,binsPhi);
601 fHistList->Add(fPtNCrossedRowsPhi);
603 fPtNCrossedRowsNClusFPhi = new TH3F("fPtNCrossedRowsNClusFPhi","fPtNCrossedRowsNClusFPhi",fgkNPtBins,binsPt,fgkNCrossedRowsNClusFBins,binsNCrossedRowsNClusF,fgkNPhiBins,binsPhi);
604 fHistList->Add(fPtNCrossedRowsNClusFPhi);
606 fPtNCrRNCrRNClusF = new TH3F("fPtNCrRNCrRNClusF","fPtNCrRNCrRNClusF",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC,fgkNCrossedRowsNClusFBins,binsNCrossedRowsNClusF);
607 fHistList->Add(fPtNCrRNCrRNClusF);
609 fPtNCrossedRowsFit = new TH2F("fPtNCrossedRowsFit","fPtNCrossedRowsFit",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC);
610 fHistList->Add(fPtNCrossedRowsFit);
612 fPtNCrossedRowsFitPhi = new TH3F("fPtNCrossedRowsFitPhi","fPtNCrossedRowsFitPhi",fgkNPtBins,binsPt,fgkNNClustersTPCBins,binsNClustersTPC,fgkNPhiBins,binsPhi);
613 fHistList->Add(fPtNCrossedRowsFitPhi);
615 fPtNCrossedRowsNClusFFitPhi = new TH3F("fPtNCrossedRowsNClusFFitPhi","fPtNCrossedRowsNClusFFitPhi",fgkNPtBins,binsPt,fgkNCrossedRowsNClusFBins,binsNCrossedRowsNClusF,fgkNPhiBins,binsPhi);
616 fHistList->Add(fPtNCrossedRowsNClusFFitPhi);
618 fNCrossedRowsNCrossedRowsFit = new TH2F("fNCrossedRowsNCrossedRowsFit","fNCrossedRowsNCrossedRowsFit",fgkNNClustersTPCBins,binsNClustersTPC,fgkNNClustersTPCBins,binsNClustersTPC);
619 fHistList->Add(fNCrossedRowsNCrossedRowsFit);
621 fNClustersNCrossedRows = new TH2F("fNClustersNCrossedRows","fNClustersNCrossedRows",fgkNNClustersTPCBins,binsNClustersTPC,fgkNNClustersTPCBins,binsNClustersTPC);
622 fHistList->Add(fNClustersNCrossedRows);
624 fNClustersNCrossedRowsFit = new TH2F("fNClustersNCrossedRowsFit","fNClustersNCrossedRowsFit",fgkNNClustersTPCBins,binsNClustersTPC,fgkNNClustersTPCBins,binsNClustersTPC);
625 fHistList->Add(fNClustersNCrossedRowsFit);
627 fPtRelUncertainty1PtNCrossedRows = new TH3F("fPtRelUncertainty1PtNCrossedRows","fPtRelUncertainty1PtNCrossedRows",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNNClustersTPCBins,binsNClustersTPC);
628 fHistList->Add(fPtRelUncertainty1PtNCrossedRows);
630 fPtRelUncertainty1PtNCrossedRowsFit = new TH3F("fPtRelUncertainty1PtNCrossedRowsFit","fPtRelUncertainty1PtNCrossedRowsFit",fgkNPtBins,binsPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty,fgkNNClustersTPCBins,binsNClustersTPC);
631 fHistList->Add(fPtRelUncertainty1PtNCrossedRowsFit);
633 fPtChi2Gold = new TH2F("fPtChi2Gold","fPtChi2Gold",fgkNPtBins,binsPt,fgkNChi2CBins,binsChi2C);
634 fHistList->Add(fPtChi2Gold);
636 fPtChi2GGC = new TH2F("fPtChi2GGC","fPtChi2GGC",fgkNPtBins,binsPt,fgkNChi2CBins,binsChi2C);
637 fHistList->Add(fPtChi2GGC);
639 fPtChi2GoldPhi = new TH3F("fPtChi2GoldPhi","fPtChi2GoldPhi",fgkNPtBins,binsPt,fgkNChi2CBins,binsChi2C,fgkNPhiBins,binsPhi);
640 fHistList->Add(fPtChi2GoldPhi);
642 fPtChi2GGCPhi = new TH3F("fPtChi2GGCPhi","fPtChi2GGCPhi",fgkNPtBins,binsPt,fgkNChi2CBins,binsChi2C,fgkNPhiBins,binsPhi);
643 fHistList->Add(fPtChi2GGCPhi);
645 fChi2GoldChi2GGC = new TH2F("fChi2GoldChi2GGC","fChi2GoldChi2GGC;#chi^{2}_{gold};#chi^{2}_{ggc}",fgkNChi2CBins,binsChi2C,fgkNChi2CBins,binsChi2C);
646 fHistList->Add(fChi2GoldChi2GGC);
649 fPtSigmaY2 = new TH2F("fPtSigmaY2","fPtSigmaY2",fgkN1PtBins,bins1Pt,fgkNSigmaY2Bins,binsSigmaY2);
650 fHistList->Add(fPtSigmaY2);
652 fPtSigmaZ2 = new TH2F("fPtSigmaZ2","fPtSigmaZ2",fgkN1PtBins,bins1Pt,fgkNSigmaZ2Bins,binsSigmaZ2);
653 fHistList->Add(fPtSigmaZ2);
655 fPtSigmaSnp2 = new TH2F("fPtSigmaSnp2","fPtSigmaSnp2",fgkN1PtBins,bins1Pt,fgkNSigmaSnp2Bins,binsSigmaSnp2);
656 fHistList->Add(fPtSigmaSnp2);
658 fPtSigmaTgl2 = new TH2F("fPtSigmaTgl2","fPtSigmaTgl2",fgkN1PtBins,bins1Pt,fgkNSigmaTgl2Bins,binsSigmaTgl2);
659 fHistList->Add(fPtSigmaTgl2);
661 fPtSigma1Pt2 = new TH2F("fPtSigma1Pt2","fPtSigma1Pt2",fgkN1PtBins,bins1Pt,fgkNSigma1Pt2Bins,binsSigma1Pt2);
662 fHistList->Add(fPtSigma1Pt2);
664 fProfPtSigmaY2 = new TProfile("fProfPtSigmaY2","fProfPtSigmaY2",fgkN1PtBins,bins1Pt);
665 fHistList->Add(fProfPtSigmaY2);
667 fProfPtSigmaZ2 = new TProfile("fProfPtSigmaZ2","fProfPtSigmaZ2",fgkN1PtBins,bins1Pt);
668 fHistList->Add(fProfPtSigmaZ2);
670 fProfPtSigmaSnp2 = new TProfile("fProfPtSigmaSnp2","fProfPtSigmaSnp2",fgkN1PtBins,bins1Pt);
671 fHistList->Add(fProfPtSigmaSnp2);
673 fProfPtSigmaTgl2 = new TProfile("fProfPtSigmaTgl2","fProfPtSigmaTgl2",fgkN1PtBins,bins1Pt);
674 fHistList->Add(fProfPtSigmaTgl2);
676 fProfPtSigma1Pt2 = new TProfile("fProfPtSigma1Pt2","fProfPtSigma1Pt2",fgkN1PtBins,bins1Pt);
677 fHistList->Add(fProfPtSigma1Pt2);
679 fProfPtSigma1Pt = new TProfile("fProfPtSigma1Pt","fProfPtSigma1Pt;p_{T};#sigma(1/p_{T})",fgkNPtBins,binsPt);
680 fHistList->Add(fProfPtSigma1Pt);
682 fProfPtPtSigma1Pt = new TProfile("fProfPtPtSigma1Pt","fProfPtPtSigma1Pt;p_{T};p_{T}#sigma(1/p_{T})",fgkNPtBins,binsPt);
683 fHistList->Add(fProfPtPtSigma1Pt);
685 TH1::AddDirectory(oldStatus);
687 PostData(1, fHistList);
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;
711 //________________________________________________________________________
712 Bool_t AliPWG4HighPtTrackQA::SelectEvent() {
714 // Decide if event should be selected for analysis
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
726 Bool_t selectEvent = kTRUE;
728 //fEvent object available?
730 AliDebug(2,Form("ERROR: fInputEvent not available\n"));
731 fNEventReject->Fill("noAliVEvent",1);
732 selectEvent = kFALSE;
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;
743 //Check if vertex is reconstructed
744 if(fDataType==kESD&&dynamic_cast<AliESDEvent*>(fEvent)) {
745 fVtx = ((AliESDEvent*)fEvent)->GetPrimaryVertexTracks();
747 if (!fVtx || !fVtx->GetStatus())
748 fVtx = ((AliESDEvent*)fEvent)->GetPrimaryVertexSPD();
751 fNEventReject->Fill("noVTX",1);
752 selectEvent = kFALSE;
756 if(!fVtx->GetStatus()) {
757 fNEventReject->Fill("VtxStatus",1);
758 selectEvent = kFALSE;
763 if(fVtx->GetNContributors()<2) {
764 fNEventReject->Fill("NCont<2",1);
765 selectEvent = kFALSE;
769 //Check if z-vertex < 10 cm
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;
778 else if(fDataType==kAOD&&dynamic_cast<AliAODEvent*>(fEvent)) {
779 const AliAODVertex *vtx = ((AliAODEvent*)fEvent)->GetPrimaryVertexSPD();
781 fNEventReject->Fill("noVTX",1);
782 selectEvent = kFALSE;
787 if(vtx->GetNContributors()<2) {
788 fNEventReject->Fill("NCont<2",1);
789 selectEvent = kFALSE;
793 //Check if z-vertex < 10 cm
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;
804 //Centrality selection should only be done in case of PbPb
807 if(fCentClass!=CalculateCentrality(fEvent) && fCentClass!=10) {
808 fNEventReject->Fill("cent",1);
809 selectEvent = kFALSE;
813 if(fDataType==kESD) {
814 if(dynamic_cast<AliESDEvent*>(fEvent)->GetCentrality()) {
815 cent = dynamic_cast<AliESDEvent*>(fEvent)->GetCentrality()->GetCentralityPercentile("V0M");
818 else if(fDataType==kAOD) {
819 if(dynamic_cast<AliAODEvent*>(fEvent)->GetHeader()->GetCentrality())
820 cent = dynamic_cast<AliAODEvent*>(fEvent)->GetHeader()->GetCentrality();
823 fNEventReject->Fill("cent>90",1);
824 selectEvent = kFALSE;
827 fh1Centrality->Fill(cent);
835 //________________________________________________________________________
836 Int_t AliPWG4HighPtTrackQA::CalculateCentrality(AliVEvent *ev){
838 // Get centrality from ESD or AOD
842 return CalculateCentrality(dynamic_cast<AliESDEvent*>(ev));
843 else if(fDataType==kAOD)
844 return CalculateCentrality(dynamic_cast<AliAODEvent*>(ev));
849 //________________________________________________________________________
850 Int_t AliPWG4HighPtTrackQA::CalculateCentrality(AliESDEvent *esd){
852 // Get centrality from ESD
858 if(esd->GetCentrality()){
859 cent = esd->GetCentrality()->GetCentralityPercentile("V0M");
860 if(fDebug>3) printf("centrality: %f\n",cent);
864 return GetCentralityClass(cent);
868 //________________________________________________________________________
869 Int_t AliPWG4HighPtTrackQA::CalculateCentrality(const AliAODEvent *aod){
871 // Get centrality from AOD
875 Float_t cent = aod->GetHeader()->GetCentrality();
876 if(fDebug>3) printf("centrality: %f\n",cent);
878 return GetCentralityClass(cent);
882 //________________________________________________________________________
883 Int_t AliPWG4HighPtTrackQA::GetCentralityClass(Float_t cent) const {
885 // Get centrality class
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;
897 //________________________________________________________________________
898 void AliPWG4HighPtTrackQA::UserExec(Option_t *) {
900 // Called for each event
901 AliDebug(2,Form(">> AliPWG4HighPtTrackQA::UserExec \n"));
903 fEvent = InputEvent();
904 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
906 // All events without selection
907 fNEventAll->Fill(0.);
911 PostData(1, fHistList);
916 //Need to keep track of selected events
917 fNEventSel->Fill(0.);
919 fVariables = new TArrayF(fNVariables);
921 if(fDataType==kESD) DoAnalysisESD();
922 if(fDataType==kAOD) DoAnalysisAOD();
924 //Delete old fVariables
925 if(fVariables) delete fVariables;
928 PostData(1, fHistList);
932 //________________________________________________________________________
933 void AliPWG4HighPtTrackQA::DoAnalysisESD() {
935 // Run analysis on ESD
939 PostData(1, fHistList);
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
947 AliMCEventHandler *eventHandlerMC = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
948 if (eventHandlerMC) {
950 if(eventHandlerMC->MCEvent()){
951 AliGenPythiaEventHeader* pythiaGenHeader = GetPythiaEventHeader(eventHandlerMC->MCEvent());
953 nTrials = pythiaGenHeader->Trials();
954 ptHard = pythiaGenHeader->GetPtHard();
956 fh1PtHard->Fill(ptHard);
957 fh1PtHardTrials->Fill(ptHard,nTrials);
959 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
964 Int_t nTracks = fESD->GetNumberOfTracks();
965 AliDebug(2,Form("nTracks ESD%d", nTracks));
968 Variables to be put in fVariables
979 10: chi2PerClusterTPC
981 12: (#crossed rows)/(#findable clusters)
987 18: NClustersTPCIter1
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
996 for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
997 fh1NTracksAll->Fill(0.);
999 //Get track for analysis
1000 AliESDtrack *track = 0x0;
1001 AliESDtrack *esdtrack = fESD->GetTrack(iTrack);
1003 fh1NTracksReject->Fill("noESDtrack",1);
1006 AliESDtrack *origtrack = new AliESDtrack(*esdtrack);
1011 if (!(fTrackCuts->AcceptTrack(esdtrack))) {
1012 fh1NTracksReject->Fill("trackCuts",1);
1013 if(origtrack) delete origtrack;
1019 track = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID());
1020 else if(fTrackType==2 || fTrackType==4) {
1021 track = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(fESD),esdtrack->GetID());
1023 fh1NTracksReject->Fill("noTPConly",1);
1024 if(origtrack) delete origtrack;
1027 AliExternalTrackParam exParam;
1028 Bool_t relate = track->RelateToVertexTPC(fVtx,fESD->GetMagneticField(),kVeryBig,&exParam);
1030 fh1NTracksReject->Fill("relate",1);
1031 if(track) delete track;
1032 if(origtrack) delete origtrack;
1035 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
1037 else if(fTrackType==5 || fTrackType==6) {
1038 if(fTrackCuts->AcceptTrack(esdtrack)) {
1039 if(origtrack) delete origtrack;
1043 if( !(fTrackCutsITSLoose->AcceptTrack(esdtrack)) && fTrackCutsTPConly->AcceptTrack(esdtrack) ) {
1046 //use TPConly constrained track
1047 track = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID());
1049 fh1NTracksReject->Fill("noTPConly",1);
1050 if(origtrack) delete origtrack;
1053 AliExternalTrackParam exParam;
1054 Bool_t relate = track->RelateToVertexTPC(fVtx,fESD->GetMagneticField(),kVeryBig,&exParam);
1056 fh1NTracksReject->Fill("relate",1);
1057 if(track) delete track;
1058 if(origtrack) delete origtrack;
1061 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
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());
1072 else if(fTrackType==7) {
1073 //use global constrained track
1074 track = new AliESDtrack(*esdtrack);
1080 if(origtrack) delete origtrack;
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;
1094 fPtAll->Fill(track->Pt());
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;
1101 if(origtrack) delete origtrack;
1106 if(fTrackCutsITSLoose ) {
1107 if(fTrackCutsITSLoose->AcceptTrack(track) ) {
1108 if(track) delete track;
1109 if(origtrack) delete origtrack;
1114 if(esdtrack->GetConstrainedParam())
1115 track->Set(esdtrack->GetConstrainedParam()->GetX(),esdtrack->GetConstrainedParam()->GetAlpha(),esdtrack->GetConstrainedParam()->GetParameter(),esdtrack->GetConstrainedParam()->GetCovariance());
1119 if(fTrackType==1 || fTrackType==2 || fTrackType==4 || fTrackType==5 || fTrackType==6 || fTrackType==7) {
1120 if(track) delete track;
1122 if(origtrack) delete origtrack;
1126 fh1NTracksSel->Fill(0.);
1128 fVariables->Reset(0.);
1130 fVariables->SetAt(track->Pt(),0);
1131 fVariables->SetAt(track->Phi(),1);
1132 fVariables->SetAt(track->Eta(),2);
1137 if(fTrackType==1 || fTrackType==2 || fTrackType==4) {
1138 track->GetImpactParametersTPC(dca2D,dcaz); //TPConly
1141 track->GetImpactParameters(dca2D,dcaz); //Global
1143 fVariables->SetAt(dca2D,3);
1144 fVariables->SetAt(dcaz,4);
1146 fVariables->SetAt((float)track->GetTPCNcls(),5);
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))
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.
1162 fVariables->SetAt(GetTrackLengthTPC(track),9);
1164 if(fVariables->At(5)>0.) fVariables->SetAt(track->GetTPCchi2()/fVariables->At(5),10);
1166 //fVariables->SetAt(track->GetTPCClusterInfo(2,1),11); //#crossed rows
1167 fVariables->SetAt(track->GetTPCCrossedRows(),11); //#crossed rows
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);
1178 fVariables->SetAt(track->GetTPCNclsIter1(),18);
1179 fVariables->SetAt(track->GetTPCchi2Iter1(),19);
1181 fVariables->SetAt(track->GetTPCnclsS(),20);
1183 Float_t chi2Gold = origtrack->GetChi2TPCConstrainedVsGlobal(fVtx);//GetGoldenChi2(origtrack);
1184 Float_t chi2GGC = GetGGCChi2(origtrack);
1186 fVariables->SetAt(chi2Gold,21);
1187 fVariables->SetAt(chi2GGC,22);
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);
1196 // int mult = fTrackCuts->CountAcceptedTracks(fESD);
1198 if(fTrackType==1 || fTrackType==2 || fTrackType==4 || fTrackType==5 || fTrackType==6 || fTrackType==7) {
1199 if(track) delete track;
1201 if(origtrack) delete origtrack;
1207 //________________________________________________________________________
1208 void AliPWG4HighPtTrackQA::DoAnalysisAOD() {
1210 // Do QA on AOD input
1212 AliAODEvent *aod = dynamic_cast<AliAODEvent*>(fEvent);
1214 AliExternalTrackParam exParam;
1215 for (Int_t iTrack = 0; iTrack < fEvent->GetNumberOfTracks(); iTrack++) {
1217 AliAODTrack *aodtrack = aod->GetTrack(iTrack);
1218 if( !aodtrack->TestFilterMask(fFilterMask) ) {
1219 fh1NTracksReject->Fill("noHybridTrack",1);
1223 fVariables->Reset(0.);
1225 fVariables->SetAt(aodtrack->Pt(),0);
1226 fVariables->SetAt(aodtrack->Phi(),1);
1227 fVariables->SetAt(aodtrack->Eta(),2);
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);
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);
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);
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);
1263 fVariables->SetAt(0.,18); //NClustersTPCIter1
1264 fVariables->SetAt(0.,19); //Chi2TPCIter1
1266 TBits sharedClusterMap = aodtrack->GetTPCSharedMap();
1267 fVariables->SetAt(sharedClusterMap.CountBits(),20);
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
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
1277 fPtAll->Fill(fVariables->At(0));
1285 //________________________________________________________________________
1286 void AliPWG4HighPtTrackQA::FillHistograms() {
1288 // Fill all QA histograms
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));
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));
1306 if(fVariables->At(18)>0.)
1307 fPtChi2PerClusterTPCIter1->Fill(fVariables->At(0),fVariables->At(19)/fVariables->At(18));
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));
1316 fPtRelUncertainty1PtITSClusterMap->Fill(fVariables->At(0),fVariables->At(0)*TMath::Sqrt(fVariables->At(17)),(int)fITSClusterMap);
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));
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)));
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)));
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));
1344 fPtChi2Gold->Fill(fVariables->At(0),fVariables->At(21));
1345 fPtChi2GGC->Fill(fVariables->At(0),fVariables->At(22));
1347 fPtChi2GoldPhi->Fill(fVariables->At(0),fVariables->At(21),fVariables->At(1));
1348 fPtChi2GGCPhi->Fill(fVariables->At(0),fVariables->At(22),fVariables->At(1));
1350 fChi2GoldChi2GGC->Fill(fVariables->At(21),fVariables->At(22));
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));
1357 fNClustersNCrossedRows->Fill(fVariables->At(5),fVariables->At(11));
1358 fNClustersNCrossedRowsFit->Fill(fVariables->At(5),fVariables->At(23));
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));
1365 //________________________________________________________________________
1366 Bool_t AliPWG4HighPtTrackQA::PythiaInfoFromFile(const char* currFile,Float_t &fXsec,Float_t &fTrials){
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
1373 TString file(currFile);
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,"");
1383 // not an archive take the basename....
1384 file.ReplaceAll(gSystem->BaseName(file.Data()),"");
1386 // Printf("%s",file.Data());
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
1391 // next trial fetch the histgram file
1392 fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
1394 // not a severe condition but inciate that we have no information
1398 // find the tlist we want to be independtent of the name so use the Tkey
1399 TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
1404 TList *list = dynamic_cast<TList*>(key->ReadObj());
1409 fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
1410 fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
1413 } // no tree pyxsec.root
1415 TTree *xtree = (TTree*)fxsec->Get("Xsection");
1421 Double_t xsection = 0;
1422 xtree->SetBranchAddress("xsection",&xsection);
1423 xtree->SetBranchAddress("ntrials",&ntrials);
1431 //________________________________________________________________________
1432 Bool_t AliPWG4HighPtTrackQA::Notify()
1435 // Implemented Notify() to read the cross sections
1436 // and number of trials from pyxsec.root
1437 // Copied from AliAnalysisTaskJetSpectrum2
1440 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
1441 Float_t xsection = 0;
1442 Float_t ftrials = 1;
1446 TFile *curfile = tree->GetCurrentFile();
1448 Error("Notify","No current file");
1451 if(!fh1Xsec||!fh1Trials){
1452 // Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
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;
1464 //________________________________________________________________________
1465 AliGenPythiaEventHeader* AliPWG4HighPtTrackQA::GetPythiaEventHeader(const AliMCEvent *mcEvent){
1467 if(!mcEvent)return 0;
1468 AliGenEventHeader* genHeader = mcEvent->GenEventHeader();
1469 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
1470 if(!pythiaGenHeader){
1472 AliGenCocktailEventHeader* genCocktailHeader = dynamic_cast<AliGenCocktailEventHeader*>(genHeader);
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__));
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)
1485 if(!pythiaGenHeader){
1486 AliWarningGeneral(Form(" %s:%d",(char*)__FILE__,__LINE__),"Pythia event header not found");
1490 return pythiaGenHeader;
1494 //_______________________________________________________________________
1495 Float_t AliPWG4HighPtTrackQA::GetTPCClusterInfo(const AliAODTrack *tr,Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1, Bool_t useFitMap) const
1497 //MV: copied from AliESDtrack since method is not available in AliAODTrack
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
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.
1511 TBits fTPCClusterMap = 0;
1513 fTPCClusterMap = tr->GetTPCFitMap();
1515 fTPCClusterMap = tr->GetTPCClusterMap();
1517 if (type==2) return fTPCClusterMap.CountBits();
1521 Int_t last=-nNeighbours;
1523 for (Int_t i=row0; i<row1; ++i){
1524 //look to current row
1525 if (fTPCClusterMap[i]) {
1531 //look to nNeighbours before
1532 if ((i-last)<=nNeighbours) {
1536 //look to nNeighbours after
1537 for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
1538 if (fTPCClusterMap[j]){
1544 if (type==1) return findable;
1549 fraction=(Float_t)found/(Float_t)findable;
1554 return 0; // undefined type - default value
1557 //_______________________________________________________________________
1558 Float_t AliPWG4HighPtTrackQA::GetTPCClusterInfoFitMap(const AliESDtrack *tr,Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1) const
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
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.
1572 TBits fTPCFitMap = tr->GetTPCFitMap();
1573 if (type==2) return fTPCFitMap.CountBits();
1577 Int_t last=-nNeighbours;
1579 for (Int_t i=row0; i<row1; ++i){
1580 //look to current row
1581 if (fTPCFitMap[i]) {
1587 //look to nNeighbours before
1588 if ((i-last)<=nNeighbours) {
1592 //look to nNeighbours after
1593 for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
1600 if (type==1) return findable;
1605 fraction=(Float_t)found/(Float_t)findable;
1610 return 0; // undefined type - default value
1613 //_______________________________________________________________________
1614 Int_t AliPWG4HighPtTrackQA::GetTrackLengthTPC(const AliESDtrack *track) const {
1616 // returns distance between 1st and last hit in TPC
1617 // distance given in number of padrows
1620 TBits fTPCClusterMap = track->GetTPCClusterMap();
1624 for(int i=0; i<=159; i++) {
1625 if(fTPCClusterMap[i]>0) firstHit = i;
1627 for(int i=159; i>=0; i--) {
1628 if(fTPCClusterMap[i]>0) lastHit = i;
1631 Int_t trackLength = lastHit - firstHit;
1636 //_______________________________________________________________________
1637 Int_t AliPWG4HighPtTrackQA::GetTrackLengthTPC(const AliAODTrack *track) const {
1639 // returns distance between 1st and last hit in TPC
1640 // distance given in number of padrows
1643 TBits fTPCClusterMap = track->GetTPCClusterMap();
1647 for(int i=0; i<=159; i++) {
1648 if(fTPCClusterMap[i]>0) firstHit = i;
1650 for(int i=159; i>=0; i--) {
1651 if(fTPCClusterMap[i]>0) lastHit = i;
1654 Int_t trackLength = lastHit - firstHit;
1659 //_______________________________________________________________________
1660 Float_t AliPWG4HighPtTrackQA::GetGoldenChi2(AliESDtrack *origtrack) {
1662 // Return chi2 between global and TPC constrained track
1663 // track should be the global unconstrained track
1666 Float_t chi2Gold = 0.;
1668 AliESDtrack *tpcTrack = 0x0;
1669 tpcTrack = AliESDtrackCuts::GetTPCOnlyTrack(fESD,origtrack->GetID());
1671 AliExternalTrackParam exParam;
1672 Bool_t relate = tpcTrack->RelateToVertexTPC(fVtx,fESD->GetMagneticField(),kVeryBig,&exParam);
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);
1678 tpcTrack->Propagate(origtrack->GetAlpha(), origtrack->GetX(), fESD->GetMagneticField());
1679 chi2Gold = (Float_t)origtrack->GetPredictedChi2(tpcTrack);
1682 if(tpcTrack) delete tpcTrack;
1688 //_______________________________________________________________________
1689 Float_t AliPWG4HighPtTrackQA::GetGGCChi2(AliESDtrack *origtrack) {
1691 // Return chi2 between global and global constrained track
1692 // track should be the global unconstrained track
1695 Float_t chi2GGC = 0.;
1697 AliESDtrack *esdtrackC = new AliESDtrack(*origtrack);
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);
1710 //________________________________________________________________________
1711 void AliPWG4HighPtTrackQA::Terminate(Option_t *)
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.