]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/AliPWG4HighPtQAMC.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Add sigma2 jet shape and fill constituent info. for subtracted jets
[u/mrichter/AliRoot.git] / PWGJE / AliPWG4HighPtQAMC.cxx
CommitLineData
fdceab34 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 compares the global reconstruction with the MC information
18// Momentum resolution is stored as function of track cuts and pt.
19// Output: Histograms for different set of cuts
20//-----------------------------------------------------------------------
21// Author : Marta Verweij - UU
22//-----------------------------------------------------------------------
23
df943115 24#ifndef ALIPWG4HighPtQAMC_CXX
25#define ALIPWG4HighPtQAMC_CXX
26
fdceab34 27#include "AliPWG4HighPtQAMC.h"
28
29#include "TH1.h"
30#include "TH2.h"
31#include "TH3.h"
b4691ee7 32#include "TProfile.h"
fdceab34 33#include "TList.h"
b4691ee7 34#include "TFile.h"
fdceab34 35#include "TChain.h"
36#include "TH3F.h"
b4691ee7 37#include "TKey.h"
38#include "TSystem.h"
39
40#include "AliAnalysisTask.h"
fdceab34 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"
df943115 49#include "AliLog.h"
b4691ee7 50#include "AliGenPythiaEventHeader.h"
51#include "AliGenCocktailEventHeader.h"
52//#include "AliAnalysisHelperJetTasks.h"
fdceab34 53
54using namespace std; //required for resolving the 'cout' symbol
55
56ClassImp(AliPWG4HighPtQAMC)
57
b4691ee7 58AliPWG4HighPtQAMC::AliPWG4HighPtQAMC()
59: AliAnalysisTask("AliPWG4HighPtQAMC", ""),
fdceab34 60 fESD(0),
b1cd0099 61 fMC(0),
cce400da 62 fStack(0),
36c36a0c 63 fVtx(0x0),
fdceab34 64 fTrackCuts(0),
327d12da 65 fTrackCutsReject(),
71e77a79 66 fTrackType(0),
0f76d8ae 67 fSigmaConstrainedMax(1e6),
b1cd0099 68 fPtMax(100.),
b4691ee7 69 fAvgTrials(1),
8f0faa80 70 fNEventAll(0),
71 fNEventSel(0),
cb76764e 72 fNEventReject(0),
b4691ee7 73 fh1Xsec(0),
74 fh1Trials(0),
75 fh1PtHard(0),
76 fh1PtHardTrials(0),
df943115 77 fPtAll(0),
78 fPtSel(0),
36c36a0c 79 fPtSelFakes(0),
80 fNPointTPCFakes(0),
81 fPtSelLargeLabel(0),
82 fMultRec(0),
83 fNPointTPCMultRec(0),
84 fDeltaPtMultRec(0),
57bee263 85 fPtAllvsPtMC(0),
a337a5a9 86 fPtAllminPtMCvsPtMC(0),
fdceab34 87 fPtAllminPtMCvsPtAll(0),
57bee263 88 fPtAllvsPtMCvsMult(0),
89 fPtAllminPtMCvsPtAllvsMult(0),
fdceab34 90 fPtAllminPtMCvsPtAllNPointTPC(0),
5a0bd31f 91 fPtAllminPtMCvsPtAllNPointTPCIter1(0),
92 fPtAllminPtMCvsPtAllChi2TPC(0),
93 fPtAllminPtMCvsPtAllChi2TPCIter1(0),
fdceab34 94 fPtAllminPtMCvsPtAllDCAR(0),
95 fPtAllminPtMCvsPtAllDCAZ(0),
96 fPtAllminPtMCvsPtAllPhi(0),
97 fPtAllminPtMCvsPtAllNPointITS(0),
98 fPtAllminPtMCvsPtAllNSigmaToVertex(0),
99 fPtAllminPtMCvsPtAllChi2C(0),
100 fPtAllminPtMCvsPtAllRel1PtUncertainty(0),
101 fPtAllMC(0),
102 fPtSelMC(0),
5a0bd31f 103 fHistList(0)
fdceab34 104{
df943115 105
b1041e3b 106 fPtBinEdges[0][0] = 10.;
107 fPtBinEdges[0][1] = 1.;
108 fPtBinEdges[1][0] = 20.;
109 fPtBinEdges[1][1] = 2.;
d889ce29 110 fPtBinEdges[2][0] = 100.;
111 fPtBinEdges[2][1] = 10.;
b1041e3b 112
fdceab34 113}
114//________________________________________________________________________
115AliPWG4HighPtQAMC::AliPWG4HighPtQAMC(const char *name):
b4691ee7 116 AliAnalysisTask(name,""),
fdceab34 117 fESD(0),
b1cd0099 118 fMC(0),
cce400da 119 fStack(0),
36c36a0c 120 fVtx(0x0),
fdceab34 121 fTrackCuts(),
327d12da 122 fTrackCutsReject(),
71e77a79 123 fTrackType(0),
0f76d8ae 124 fSigmaConstrainedMax(1e6),
b1cd0099 125 fPtMax(100.),
b4691ee7 126 fAvgTrials(1),
8f0faa80 127 fNEventAll(0),
128 fNEventSel(0),
cb76764e 129 fNEventReject(0),
b4691ee7 130 fh1Xsec(0),
131 fh1Trials(0),
132 fh1PtHard(0),
133 fh1PtHardTrials(0),
df943115 134 fPtAll(0),
135 fPtSel(0),
36c36a0c 136 fPtSelFakes(0),
137 fNPointTPCFakes(0),
138 fPtSelLargeLabel(0),
139 fMultRec(0),
140 fNPointTPCMultRec(0),
141 fDeltaPtMultRec(0),
57bee263 142 fPtAllvsPtMC(0),
a337a5a9 143 fPtAllminPtMCvsPtMC(0),
fdceab34 144 fPtAllminPtMCvsPtAll(0),
57bee263 145 fPtAllvsPtMCvsMult(0),
146 fPtAllminPtMCvsPtAllvsMult(0),
fdceab34 147 fPtAllminPtMCvsPtAllNPointTPC(0),
5a0bd31f 148 fPtAllminPtMCvsPtAllNPointTPCIter1(0),
149 fPtAllminPtMCvsPtAllChi2TPC(0),
150 fPtAllminPtMCvsPtAllChi2TPCIter1(0),
fdceab34 151 fPtAllminPtMCvsPtAllDCAR(0),
152 fPtAllminPtMCvsPtAllDCAZ(0),
153 fPtAllminPtMCvsPtAllPhi(0),
154 fPtAllminPtMCvsPtAllNPointITS(0),
155 fPtAllminPtMCvsPtAllNSigmaToVertex(0),
156 fPtAllminPtMCvsPtAllChi2C(0),
157 fPtAllminPtMCvsPtAllRel1PtUncertainty(0),
158 fPtAllMC(0),
159 fPtSelMC(0),
5a0bd31f 160 fHistList(0)
fdceab34 161{
162 //
163 // Constructor. Initialization of Inputs and Outputs
164 //
f51451be 165 AliDebug(2,Form("AliPWG4HighPtQAMC Calling Constructor"));
b1041e3b 166
167 fPtBinEdges[0][0] = 10.;
168 fPtBinEdges[0][1] = 1.;
169 fPtBinEdges[1][0] = 20.;
170 fPtBinEdges[1][1] = 2.;
d889ce29 171 fPtBinEdges[2][0] = 100.;
172 fPtBinEdges[2][1] = 10.;
b1041e3b 173
fdceab34 174 // Input slot #0 works with a TChain ESD
175 DefineInput(0, TChain::Class());
b4691ee7 176 // Output slot #0, #1 write into a TList
fdceab34 177 DefineOutput(0, TList::Class());
df943115 178}
fdceab34 179
b1041e3b 180//________________________________________________________________________
181void AliPWG4HighPtQAMC::SetPtBinEdges(Int_t region, Double_t ptmax, Double_t ptBinWidth) {
182
183 if(region<3) {
184 fPtBinEdges[region][0] = ptmax;
185 fPtBinEdges[region][1] = ptBinWidth;
186 }
187 else {
188 AliError("Only 3 regions alowed. Use region 0/1/2\n");
189 return;
190 }
191
192}
193
fdceab34 194//________________________________________________________________________
195void AliPWG4HighPtQAMC::ConnectInputData(Option_t *)
196{
b1cd0099 197 // Connect ESD and MC event handler here
fdceab34 198 // Called once
df943115 199 AliDebug(2,Form(">> AliPWG4HighPtSpectra::ConnectInputData \n"));
fdceab34 200 TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
201 if (!tree) {
b1cd0099 202 AliDebug(2,Form( "ERROR: Could not read chain from input slot 0 \n"));
203 return;
fdceab34 204 }
b1cd0099 205
206 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
207
208 if (!esdH) {
209 AliDebug(2,Form("ERROR: Could not get ESDInputHandler"));
210 return;
211 } else
212 fESD = esdH->GetEvent();
213
cce400da 214 AliMCEventHandler *eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
215 if (!eventHandler) {
216 AliDebug(2,Form( "ERROR: Could not retrieve MC event handler \n"));
217 }
b1cd0099 218 else
219 fMC = eventHandler->MCEvent();
220
fdceab34 221}
222
b4691ee7 223
fdceab34 224//________________________________________________________________________
225void AliPWG4HighPtQAMC::CreateOutputObjects() {
226 //Create output objects
df943115 227 AliDebug(2,Form(">> AliPWG4HighPtQATPConly::CreateOutputObjects \n"));
228
229 Bool_t oldStatus = TH1::AddDirectoryStatus();
230 TH1::AddDirectory(kFALSE);
231
fdceab34 232 OpenFile(0);
233 fHistList = new TList();
b4691ee7 234 fHistList->SetOwner(kTRUE);
fdceab34 235
b1041e3b 236 Float_t fgkPtMin = 0.;
237
238 //fPtBinEdges[region][0] = ptmax of region ; fPtBinEdges[region][1] = binWidth of region
239 const Float_t ptmin1 = fgkPtMin;
240 const Float_t ptmax1 = fPtBinEdges[0][0];
241 const Float_t ptmin2 = ptmax1 ;
242 const Float_t ptmax2 = fPtBinEdges[1][0];
243 const Float_t ptmin3 = ptmax2 ;
244 const Float_t ptmax3 = fPtBinEdges[2][0];//fgkPtMax;
245 const Int_t nbin11 = (int)((ptmax1-ptmin1)/fPtBinEdges[0][1]);
246 const Int_t nbin12 = (int)((ptmax2-ptmin2)/fPtBinEdges[1][1])+nbin11;
247 const Int_t nbin13 = (int)((ptmax3-ptmin3)/fPtBinEdges[2][1])+nbin12;
248 Int_t fgkNPtBins=nbin13;
249 //Create array with low edges of each bin
250 Double_t *binsPt=new Double_t[fgkNPtBins+1];
251 for(Int_t i=0; i<=fgkNPtBins; i++) {
252 if(i<=nbin11) binsPt[i]=(Double_t)ptmin1 + (ptmax1-ptmin1)/nbin11*(Double_t)i ;
253 if(i<=nbin12 && i>nbin11) binsPt[i]=(Double_t)ptmin2 + (ptmax2-ptmin2)/(nbin12-nbin11)*((Double_t)i-(Double_t)nbin11) ;
254 if(i<=nbin13 && i>nbin12) binsPt[i]=(Double_t)ptmin3 + (ptmax3-ptmin3)/(nbin13-nbin12)*((Double_t)i-(Double_t)nbin12) ;
255 }
256
257 Int_t fgkNPhiBins = 18*6;
258 Float_t kMinPhi = 0.;
259 Float_t kMaxPhi = 2.*TMath::Pi();
260 Double_t *binsPhi = new Double_t[fgkNPhiBins+1];
261 for(Int_t i=0; i<=fgkNPhiBins; i++) binsPhi[i]=(Double_t)kMinPhi + (kMaxPhi-kMinPhi)/fgkNPhiBins*(Double_t)i ;
262
8f0faa80 263 Int_t fgkResPtBins=80;
b1041e3b 264 Float_t kMinResPt = -1.;
265 Float_t kMaxResPt = 1.;
266 Double_t *binsResPt = new Double_t[fgkResPtBins+1];
267 for(Int_t i=0; i<=fgkResPtBins; i++) binsResPt[i]=(Double_t)kMinResPt + (kMaxResPt-kMinResPt)/fgkResPtBins*(Double_t)i ;
268
269 Int_t fgkMultBins=50;
270 Float_t kMinMult = 0.;
271 Float_t kMaxMult = 3000.;
272 Double_t *binsMult = new Double_t[fgkMultBins+1];
273 for(Int_t i=0; i<=fgkMultBins; i++) binsMult[i]=(Double_t)kMinMult + (kMaxMult-kMinMult)/fgkMultBins*(Double_t)i ;
274
275 Int_t fgkNEtaBins=20;
276 Float_t fgkEtaMin = -1.;
277 Float_t fgkEtaMax = 1.;
278 Double_t *binsEta=new Double_t[fgkNEtaBins+1];
279 for(Int_t i=0; i<=fgkNEtaBins; i++) binsEta[i]=(Double_t)fgkEtaMin + (fgkEtaMax-fgkEtaMin)/fgkNEtaBins*(Double_t)i ;
280
281 Int_t fgkNNClustersTPCBins=80;
282 Float_t fgkNClustersTPCMin = 0.5;
283 Float_t fgkNClustersTPCMax = 160.5;
284 Double_t *binsNClustersTPC=new Double_t[fgkNNClustersTPCBins+1];
285 for(Int_t i=0; i<=fgkNNClustersTPCBins; i++) binsNClustersTPC[i]=(Double_t)fgkNClustersTPCMin + (fgkNClustersTPCMax-fgkNClustersTPCMin)/fgkNNClustersTPCBins*(Double_t)i ;
286
287 Int_t fgkNDCA2DBins=80;
288 Float_t fgkDCA2DMin = -0.2;
289 Float_t fgkDCA2DMax = 0.2;
290 if(fTrackType==1 || fTrackType==2 || fTrackType==4) {
291 fgkDCA2DMin = -2.;
292 fgkDCA2DMax = 2.;
293 }
294 Double_t *binsDCA2D=new Double_t[fgkNDCA2DBins+1];
295 for(Int_t i=0; i<=fgkNDCA2DBins; i++) binsDCA2D[i]=(Double_t)fgkDCA2DMin + (fgkDCA2DMax-fgkDCA2DMin)/fgkNDCA2DBins*(Double_t)i ;
296
297 Int_t fgkNDCAZBins=80;
298 Float_t fgkDCAZMin = -2.;
299 Float_t fgkDCAZMax = 2.;
300 if(fTrackType==1 || fTrackType==2 || fTrackType==4) {
301 fgkDCAZMin = -5.;
302 fgkDCAZMax = 5.;
303 }
304 Double_t *binsDCAZ=new Double_t[fgkNDCAZBins+1];
305 for(Int_t i=0; i<=fgkNDCAZBins; i++) binsDCAZ[i]=(Double_t)fgkDCAZMin + (fgkDCAZMax-fgkDCAZMin)/fgkNDCAZBins*(Double_t)i ;
306
307 Int_t fgkNNPointITSBins=9;
308 Float_t fgkNPointITSMin = -0.5;
309 Float_t fgkNPointITSMax = 8.5;
310 Double_t *binsNPointITS=new Double_t[fgkNNPointITSBins+1];
311 for(Int_t i=0; i<=fgkNNPointITSBins; i++) binsNPointITS[i]=(Double_t)fgkNPointITSMin + (fgkNPointITSMax-fgkNPointITSMin)/fgkNNPointITSBins*(Double_t)i ;
312
313 Int_t fgkNNSigmaToVertexBins=20;
314 Float_t fgkNSigmaToVertexMin = 0.;
315 Float_t fgkNSigmaToVertexMax = 8.;
316 Double_t *binsNSigmaToVertex=new Double_t[fgkNNSigmaToVertexBins+1];
317 for(Int_t i=0; i<=fgkNNSigmaToVertexBins; i++) binsNSigmaToVertex[i]=(Double_t)fgkNSigmaToVertexMin + (fgkNSigmaToVertexMax-fgkNSigmaToVertexMin)/fgkNNSigmaToVertexBins*(Double_t)i ;
318
319 Int_t fgkNChi2CBins=20;
320 Float_t fgkChi2CMin = 0.;
321 Float_t fgkChi2CMax = 100.;
322 Double_t *binsChi2C=new Double_t[fgkNChi2CBins+1];
323 for(Int_t i=0; i<=fgkNChi2CBins; i++) binsChi2C[i]=(Double_t)fgkChi2CMin + (fgkChi2CMax-fgkChi2CMin)/fgkNChi2CBins*(Double_t)i ;
324
325 Int_t fgkNRel1PtUncertaintyBins=50;
326 Float_t fgkRel1PtUncertaintyMin = 0.;
327 Float_t fgkRel1PtUncertaintyMax = 1.;
328 Double_t *binsRel1PtUncertainty=new Double_t[fgkNRel1PtUncertaintyBins+1];
329 for(Int_t i=0; i<=fgkNRel1PtUncertaintyBins; i++) binsRel1PtUncertainty[i]=(Double_t)fgkRel1PtUncertaintyMin + (fgkRel1PtUncertaintyMax-fgkRel1PtUncertaintyMin)/fgkNRel1PtUncertaintyBins*(Double_t)i ;
330
331 Float_t fgkChi2PerClusMin = 0.;
332 Float_t fgkChi2PerClusMax = 4.;
333 Int_t fgkNChi2PerClusBins = (int)(fgkChi2PerClusMax*10.);
334 Double_t *binsChi2PerClus=new Double_t[fgkNChi2PerClusBins+1];
335 for(Int_t i=0; i<=fgkNChi2PerClusBins; i++) binsChi2PerClus[i]=(Double_t)fgkChi2PerClusMin + (fgkChi2PerClusMax-fgkChi2PerClusMin)/fgkNChi2PerClusBins*(Double_t)i ;
fdceab34 336
57bee263 337
8f0faa80 338 fNEventAll = new TH1F("fNEventAll","NEventAll",1,-0.5,0.5);
339 fHistList->Add(fNEventAll);
340 fNEventSel = new TH1F("fNEventSel","NEvent Selected for analysis",1,-0.5,0.5);
341 fHistList->Add(fNEventSel);
cb76764e 342 fNEventReject = new TH1F("fNEventReject","Reason events are rejectected for analysis",20,0,20);
36c36a0c 343 //Set labels
344 fNEventReject->Fill("noESD",0);
345 fNEventReject->Fill("Trigger",0);
346 fNEventReject->Fill("noMCEvent",0);
11245558 347 fNEventReject->Fill("noStack",0);
36c36a0c 348 fNEventReject->Fill("NTracks<2",0);
349 fNEventReject->Fill("noVTX",0);
350 fNEventReject->Fill("VtxStatus",0);
351 fNEventReject->Fill("NCont<2",0);
352 fNEventReject->Fill("ZVTX>10",0);
cb76764e 353 fHistList->Add(fNEventReject);
b4691ee7 354
355 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
356 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
357 fHistList->Add(fh1Xsec);
358
359 fh1Trials = new TH1F("fh1Trials","trials root file",1,0,1);
360 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
361 fHistList->Add(fh1Trials);
362
363 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
364 fHistList->Add(fh1PtHard);
365 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
366 fHistList->Add(fh1PtHardTrials);
367
b1041e3b 368 fPtAll = new TH1F("fPtAll","PtAll",fgkNPtBins, binsPt);
fdceab34 369 fHistList->Add(fPtAll);
b1041e3b 370 fPtSel = new TH1F("fPtSel","PtSel",fgkNPtBins, binsPt);
fdceab34 371 fHistList->Add(fPtSel);
36c36a0c 372
b1041e3b 373 fPtSelFakes = new TH1F("fPtSelFakes","PtSelFakes",fgkNPtBins, binsPt);
36c36a0c 374 fHistList->Add(fPtSelFakes);
375
b1041e3b 376 fNPointTPCFakes = new TH1F("fNPointTPCFakes","fNPointTPCFakes",fgkNNClustersTPCBins, binsNClustersTPC);
36c36a0c 377 fHistList->Add(fNPointTPCFakes);
378
b1041e3b 379 fPtSelLargeLabel = new TH1F("fPtSelLargeLabel","PtSelLargeLabel",fgkNPtBins, binsPt);
36c36a0c 380 fHistList->Add(fPtSelLargeLabel);
381
b1041e3b 382 fMultRec = new TH1F("fMultRec","Multiple reconstruction of tracks",fgkMultBins, binsMult);
36c36a0c 383 fHistList->Add(fMultRec);
384
b1041e3b 385 fNPointTPCMultRec = new TH1F("fNPointTPCMultRec","Multiple reconstruction of tracks",fgkNNClustersTPCBins, binsNClustersTPC);
36c36a0c 386 fHistList->Add(fNPointTPCMultRec);
387
b1041e3b 388 fDeltaPtMultRec = new TH2F("fDeltaPtMultRec","Delta pT vs pT for multiple reconstructed tracks",50,0.,50.,40,-20.,20.);
36c36a0c 389 fHistList->Add(fDeltaPtMultRec);
390
b1041e3b 391 fPtAllvsPtMC = new TH2F("fPtAllvsPtMC","fPtAllvsPtMC;p_{T,MC};p_{T,rec}",fgkNPtBins, binsPt,fgkNPtBins, binsPt);
57bee263 392 fHistList->Add(fPtAllvsPtMC);
393
a337a5a9 394 fPtAllminPtMCvsPtMC = new TH2F("fPtAllminPtMCvsPtMC","PtAllminPtMCvsPtMC",fgkNPtBins, binsPt,fgkResPtBins,binsResPt);
395 fPtAllminPtMCvsPtMC->SetXTitle("p_{t}^{MC}");
396 fPtAllminPtMCvsPtMC->SetYTitle("(1/p_{t}^{rec}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
397 fHistList->Add(fPtAllminPtMCvsPtMC);
398
b1041e3b 399 fPtAllminPtMCvsPtAll = new TH2F("fPtAllminPtMCvsPtAll","PtAllminPtMCvsPtAll",fgkNPtBins, binsPt,fgkResPtBins,binsResPt);
a337a5a9 400 fPtAllminPtMCvsPtAll->SetXTitle("p_{t}^{rec}");
401 fPtAllminPtMCvsPtAll->SetYTitle("(1/p_{t}^{rec}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
fdceab34 402 fHistList->Add(fPtAllminPtMCvsPtAll);
57bee263 403
b1041e3b 404 fPtAllvsPtMCvsMult = new TH3F("fPtAllvsPtMCvsMult","fPtAllvsPtMCvsMult;p_{T,MC};p_{T,rec};N_{tracks}",fgkNPtBins, binsPt,fgkNPtBins, binsPt,fgkMultBins,binsMult);
57bee263 405 fHistList->Add(fPtAllvsPtMCvsMult);
406
b1041e3b 407 fPtAllminPtMCvsPtAllvsMult = new TH3F("fPtAllminPtMCvsPtAllvsMult","fPtAllminPtMCvsPtAllvsMult",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkMultBins,binsMult);
57bee263 408 fPtAllminPtMCvsPtAllvsMult->SetXTitle("p_{t}^{MC}");
409 fPtAllminPtMCvsPtAllvsMult->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
410 fPtAllminPtMCvsPtAllvsMult->SetZTitle("N_{tracks}");
411 fHistList->Add(fPtAllminPtMCvsPtAllvsMult);
fdceab34 412
b1041e3b 413 fPtAllminPtMCvsPtAllNPointTPC = new TH3F("fPtAllminPtMCvsPtAllNPointTPC","PtAllminPtMCvsPtAllNPointTPC",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkNNClustersTPCBins, binsNClustersTPC);
fdceab34 414 fPtAllminPtMCvsPtAllNPointTPC->SetXTitle("p_{t}^{MC}");
415 fPtAllminPtMCvsPtAllNPointTPC->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
5a0bd31f 416 fPtAllminPtMCvsPtAllNPointTPC->SetZTitle("N_{cls,TPC}");
fdceab34 417 fHistList->Add(fPtAllminPtMCvsPtAllNPointTPC);
418
b1041e3b 419 fPtAllminPtMCvsPtAllNPointTPCIter1 = new TH3F("fPtAllminPtMCvsPtAllNPointTPCIter1","PtAllminPtMCvsPtAllNPointTPCIter1",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkNNClustersTPCBins, binsNClustersTPC);
5a0bd31f 420 fPtAllminPtMCvsPtAllNPointTPCIter1->SetXTitle("p_{t}^{MC}");
421 fPtAllminPtMCvsPtAllNPointTPCIter1->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
422 fPtAllminPtMCvsPtAllNPointTPCIter1->SetZTitle("N_{clsIter1,TPC}");
423 fHistList->Add(fPtAllminPtMCvsPtAllNPointTPCIter1);
424
b1041e3b 425 fPtAllminPtMCvsPtAllChi2TPC = new TH3F("fPtAllminPtMCvsPtAllChi2TPC","PtAllminPtMCvsPtAllChi2TPC",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkNChi2PerClusBins, binsChi2PerClus);
5a0bd31f 426 fPtAllminPtMCvsPtAllChi2TPC->SetXTitle("p_{t}^{MC}");
427 fPtAllminPtMCvsPtAllChi2TPC->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
428 fPtAllminPtMCvsPtAllChi2TPC->SetZTitle("#chi^{2} TPC");
429 fHistList->Add(fPtAllminPtMCvsPtAllChi2TPC);
430
b1041e3b 431 fPtAllminPtMCvsPtAllChi2TPCIter1 = new TH3F("fPtAllminPtMCvsPtAllChi2TPCIter1","PtAllminPtMCvsPtAllChi2TPCIter1",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkNChi2PerClusBins, binsChi2PerClus);
5a0bd31f 432 fPtAllminPtMCvsPtAllChi2TPCIter1->SetXTitle("p_{t}^{MC}");
433 fPtAllminPtMCvsPtAllChi2TPCIter1->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
434 fPtAllminPtMCvsPtAllChi2TPCIter1->SetZTitle("#chi^{2} TPC Iter1");
435 fHistList->Add(fPtAllminPtMCvsPtAllChi2TPCIter1);
436
b1041e3b 437 fPtAllminPtMCvsPtAllDCAR = new TH3F("fPtAllminPtMCvsPtAllDCAR","PtAllminPtMCvsPtAllDCAR",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkNDCA2DBins,binsDCA2D);
fdceab34 438 fPtAllminPtMCvsPtAllDCAR->SetXTitle("p_{t}^{MC}");
439 fPtAllminPtMCvsPtAllDCAR->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
440 fPtAllminPtMCvsPtAllDCAR->SetZTitle("DCA_{R}");
441 fHistList->Add(fPtAllminPtMCvsPtAllDCAR);
442
b1041e3b 443 fPtAllminPtMCvsPtAllDCAZ = new TH3F("fPtAllminPtMCvsPtAllDCAZ","PtAllminPtMCvsPtAllDCAZ",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkNDCAZBins,binsDCAZ);
fdceab34 444 fPtAllminPtMCvsPtAllDCAZ->SetXTitle("p_{t}^{MC}");
445 fPtAllminPtMCvsPtAllDCAZ->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
446 fPtAllminPtMCvsPtAllDCAZ->SetZTitle("DCA_{Z}");
447 fHistList->Add(fPtAllminPtMCvsPtAllDCAZ);
448
b1041e3b 449 fPtAllminPtMCvsPtAllPhi = new TH3F("fPtAllminPtMCvsPtAllPhi","PtAllminPtMCvsPtAllPhi",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkNPhiBins,binsPhi);
fdceab34 450 fPtAllminPtMCvsPtAllPhi->SetXTitle("p_{t}^{MC}");
451 fPtAllminPtMCvsPtAllPhi->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
452 fPtAllminPtMCvsPtAllPhi->SetZTitle("#phi");
453 fHistList->Add(fPtAllminPtMCvsPtAllPhi);
454
b1041e3b 455 fPtAllminPtMCvsPtAllNPointITS = new TH3F("fPtAllminPtMCvsPtAllNPointITS","PtAllminPtMCvsPtAllNPointITS",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkNNPointITSBins,binsNPointITS);
fdceab34 456 fPtAllminPtMCvsPtAllNPointITS->SetXTitle("p_{t}^{MC}");
457 fPtAllminPtMCvsPtAllNPointITS->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
458 fPtAllminPtMCvsPtAllNPointITS->SetZTitle("N_{point,ITS}}");
459 fHistList->Add(fPtAllminPtMCvsPtAllNPointITS);
460
b1041e3b 461 fPtAllminPtMCvsPtAllNSigmaToVertex = new TH3F("fPtAllminPtMCvsPtAllNSigmaToVertex","PtAllminPtMCvsPtAllNSigmaToVertex",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkNNSigmaToVertexBins,binsNSigmaToVertex);
fdceab34 462 fPtAllminPtMCvsPtAllNSigmaToVertex->SetXTitle("p_{t}^{MC}");
463 fPtAllminPtMCvsPtAllNSigmaToVertex->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
464 fPtAllminPtMCvsPtAllNSigmaToVertex->SetZTitle("N#sigma to vertex");
465 fHistList->Add(fPtAllminPtMCvsPtAllNSigmaToVertex);
466
b1041e3b 467 fPtAllminPtMCvsPtAllChi2C = new TH3F("fPtAllminPtMCvsPtAllChi2C","PtAllminPtMCvsPtAllChi2C",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkNChi2CBins,binsChi2C);
fdceab34 468 fPtAllminPtMCvsPtAllChi2C->SetXTitle("p_{t}^{MC}");
469 fPtAllminPtMCvsPtAllChi2C->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
470 fPtAllminPtMCvsPtAllChi2C->SetZTitle("Constrained #chi^{2}");
471 fHistList->Add(fPtAllminPtMCvsPtAllChi2C);
472
b1041e3b 473 fPtAllminPtMCvsPtAllRel1PtUncertainty = new TH3F("fPtAllminPtMCvsPtAllRel1PtUncertainty","PtAllminPtMCvsPtAllRel1PtUncertainty",fgkNPtBins, binsPt,fgkResPtBins,binsResPt,fgkNRel1PtUncertaintyBins,binsRel1PtUncertainty);
fdceab34 474 fPtAllminPtMCvsPtAllRel1PtUncertainty->SetXTitle("p_{t}^{MC}");
475 fPtAllminPtMCvsPtAllRel1PtUncertainty->SetYTitle("(1/p_{t}^{All}-1/p_{t}^{MC})/(1/p_{t}^{MC})");
476 fPtAllminPtMCvsPtAllRel1PtUncertainty->SetZTitle("Rel1PtUncertainty");
477 fHistList->Add(fPtAllminPtMCvsPtAllRel1PtUncertainty);
478
fdceab34 479
b1041e3b 480 fPtAllMC = new TH1F("fPtAllMC","PtAll",fgkNPtBins, binsPt);
fdceab34 481 fHistList->Add(fPtAllMC);
b1041e3b 482 fPtSelMC = new TH1F("fPtSelMC","PtSel",fgkNPtBins, binsPt);
fdceab34 483 fHistList->Add(fPtSelMC);
b1041e3b 484
df943115 485 TH1::AddDirectory(oldStatus);
486
11245558 487 PostData(0, fHistList);
b1041e3b 488
489 if(binsPt) delete [] binsPt;
490 if(binsResPt) delete [] binsResPt;
491 if(binsPhi) delete [] binsPhi;
492 if(binsNClustersTPC) delete [] binsNClustersTPC;
493 if(binsDCA2D) delete [] binsDCA2D;
494 if(binsDCAZ) delete [] binsDCAZ;
495 if(binsNPointITS) delete [] binsNPointITS;
496 if(binsNSigmaToVertex) delete [] binsNSigmaToVertex;
497 if(binsChi2C) delete [] binsChi2C;
498 if(binsEta) delete [] binsEta;
499 if(binsRel1PtUncertainty) delete [] binsRel1PtUncertainty;
500 if(binsChi2PerClus) delete [] binsChi2PerClus;
501 if(binsMult) delete [] binsMult;
502
fdceab34 503}
cce400da 504
fdceab34 505//________________________________________________________________________
cce400da 506Bool_t AliPWG4HighPtQAMC::SelectEvent() {
507 //
508 // Decide if event should be selected for analysis
509 //
8f0faa80 510
cce400da 511 // Checks following requirements:
512 // - fESD available
513 // - trigger info from AliPhysicsSelection
514 // - MCevent available
515 // - number of reconstructed tracks > 1
516 // - primary vertex reconstructed
517 // - z-vertex < 10 cm
518
519 Bool_t selectEvent = kTRUE;
520
521 //fESD object available?
fdceab34 522 if (!fESD) {
b4691ee7 523 AliDebug(2,Form("ERROR: fInputEvent not available\n"));
cb76764e 524 fNEventReject->Fill("noESD",1);
cce400da 525 selectEvent = kFALSE;
526 return selectEvent;
fdceab34 527 }
528
cce400da 529 //Trigger
bfd58011 530 UInt_t isSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
531 if(!(isSelected&AliVEvent::kMB)) { //Select collison candidates
8f0faa80 532 AliDebug(2,Form(" Trigger Selection: event REJECTED ... "));
cb76764e 533 fNEventReject->Fill("Trigger",1);
cce400da 534 selectEvent = kFALSE;
535 return selectEvent;
8f0faa80 536 }
cce400da 537
538 //MCEvent available?
539 //if yes: get stack
b1cd0099 540 if(fMC) {
541 AliDebug(2,Form("MC particles: %d", fMC->GetNumberOfTracks()));
cce400da 542 fStack = fMC->Stack(); //Particles Stack
11245558 543 if(fStack) {
544 AliDebug(2,Form("MC particles stack: %d", fStack->GetNtrack()));
545 }
546 else {
547 AliDebug(2,Form("ERROR: Could not retrieve MC eventHandler"));
548 fNEventReject->Fill("noStack",1);
549 selectEvent = kFALSE;
550 return selectEvent;
551 }
b1cd0099 552 } else {
11245558 553 AliDebug(2,Form("ERROR: Could not retrieve stack"));
cb76764e 554 fNEventReject->Fill("noMCEvent",1);
cce400da 555 selectEvent = kFALSE;
556 return selectEvent;
fdceab34 557 }
558
cce400da 559 //Check if number of reconstructed tracks is larger than 1
560 if(!fESD->GetNumberOfTracks() || fESD->GetNumberOfTracks()<2) {
561 fNEventReject->Fill("NTracks<2",1);
562 selectEvent = kFALSE;
563 return selectEvent;
cb76764e 564 }
b4691ee7 565
cce400da 566 //Check if vertex is reconstructed
36c36a0c 567 if(fTrackType==1) fVtx = fESD->GetPrimaryVertexTPC();
568 else fVtx = fESD->GetPrimaryVertexSPD();
569
570 if(!fVtx) {
571 fNEventReject->Fill("noVTX",1);
572 selectEvent = kFALSE;
573 return selectEvent;
574 }
575
576 if(!fVtx->GetStatus()) {
577 fNEventReject->Fill("VtxStatus",1);
578 selectEvent = kFALSE;
579 return selectEvent;
580 }
581
fdceab34 582 // Need vertex cut
36c36a0c 583 // TString vtxName(fVtx->GetName());
584 if(fVtx->GetNContributors()<2) {
585 fNEventReject->Fill("NCont<2",1);
586 selectEvent = kFALSE;
587 return selectEvent;
8f0faa80 588 }
36c36a0c 589
cce400da 590 //Check if z-vertex < 10 cm
75946d5d 591 double primVtx[3];
36c36a0c 592 fVtx->GetXYZ(primVtx);
8f0faa80 593 if(TMath::Sqrt(primVtx[0]*primVtx[0] + primVtx[1]*primVtx[1])>1. || TMath::Abs(primVtx[2]>10.)){
cb76764e 594 fNEventReject->Fill("ZVTX>10",1);
cce400da 595 selectEvent = kFALSE;
596 return selectEvent;
8f0faa80 597 }
36c36a0c 598
599 AliDebug(2,Form("Vertex title %s, status %d, nCont %d\n",fVtx->GetTitle(), fVtx->GetStatus(), fVtx->GetNContributors()));
fdceab34 600
cce400da 601 return selectEvent;
602
603}
604
605//________________________________________________________________________
606void AliPWG4HighPtQAMC::Exec(Option_t *) {
607 // Main loop
608 // Called for each event
609 AliDebug(2,Form(">> AliPWG4HighPtQATPConly::Exec \n"));
610 // All events without selection
611 fNEventAll->Fill(0.);
612
613 if(!SelectEvent()) {
614 // Post output data
8f0faa80 615 PostData(0, fHistList);
8f0faa80 616 return;
617 }
fdceab34 618
cce400da 619 // ---- Get MC Header information (for MC productions in pThard bins) ----
620 Double_t ptHard = 0.;
621 Double_t nTrials = 1; // trials for MC trigger weight for real data
622
623 if(fMC){
624 AliGenPythiaEventHeader* pythiaGenHeader = GetPythiaEventHeader(fMC);
625 if(pythiaGenHeader){
626 nTrials = pythiaGenHeader->Trials();
627 ptHard = pythiaGenHeader->GetPtHard();
628
629 fh1PtHard->Fill(ptHard);
630 fh1PtHardTrials->Fill(ptHard,nTrials);
631
632 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
633 }
634 }
635
b4691ee7 636 //Need to keep track of selected events
637 fNEventSel->Fill(0.);
638
fdceab34 639 Int_t nTracks = fESD->GetNumberOfTracks();
b1cd0099 640 AliDebug(2,Form("nTracks ESD%d", nTracks));
df943115 641
cce400da 642 int nMCtracks = fStack->GetNtrack();
cba109a0 643
71e77a79 644 Float_t pt = 0.;
645 Float_t ptMC = 0.;
646 Float_t phi = 0.;
647 Float_t dca2D = 0.;
648 Float_t dcaZ = 0.;
649 Int_t nPointITS = 0;
650 Float_t chi2C = 0.;
651 Float_t nSigmaToVertex = 0.;
652 Float_t relUncertainty1Pt = 0.;
653
13e88539 654 int mult = fTrackCuts->CountAcceptedTracks(fESD);
655
fdceab34 656 for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
657
11cb9d82 658 AliESDtrack *track = 0;
71e77a79 659 AliESDtrack *esdtrack = fESD->GetTrack(iTrack);
660 if(!esdtrack) continue;
71e77a79 661
8da84c2c 662 if(fTrackType==1) {
71e77a79 663 track = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID());
8da84c2c 664 if(!track) continue;
665 }
36c36a0c 666 else if(fTrackType==2) {
667 track = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID());
8da84c2c 668 if(!track) continue;
669
36c36a0c 670 AliExternalTrackParam exParam;
671 Bool_t relate = track->RelateToVertexTPC(fVtx,fESD->GetMagneticField(),kVeryBig,&exParam);
672 if( !relate ) {
673 delete track;
674 continue;
675 }
676 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
677 }
d3a3f33d 678 else if(fTrackType==7) {
679 //use global constrained track
42881dab 680 track = new AliESDtrack(*esdtrack);
681 // track = esdtrack;
031aac8e 682 // track->Set(esdtrack->GetConstrainedParam()->GetX(),esdtrack->GetConstrainedParam()->GetAlpha(),esdtrack->GetConstrainedParam()->GetParameter(),esdtrack->GetConstrainedParam()->GetCovariance());
d3a3f33d 683 }
71e77a79 684 else
685 track = esdtrack;
83ff93cd 686
687 if(!track) continue;
71e77a79 688
aa3ba8d2 689 if(fTrackType==2) {
690 //Cut on chi2 of constrained fit
691 if(track->GetConstrainedChi2TPC() > fSigmaConstrainedMax*fSigmaConstrainedMax) {
692 delete track;
693 continue;
694 }
695 }
696
fdceab34 697 Int_t label = TMath::Abs(track->GetLabel());
36c36a0c 698 if(label>=nMCtracks) {
699 if (fTrackCuts->AcceptTrack(track)) { fPtSelLargeLabel->Fill(pt); }
42881dab 700 if(fTrackType==1 || fTrackType==2 || fTrackType==7) delete track;
36c36a0c 701 continue;
702 }
031aac8e 703
cce400da 704 TParticle *particle = fStack->Particle(label) ;
36c36a0c 705 if(!particle) {
42881dab 706 if(fTrackType==1 || fTrackType==2 || fTrackType==7) {
b1041e3b 707 if(track) delete track;
708 }
36c36a0c 709 continue;
710 }
fdceab34 711
71e77a79 712 ptMC = particle->Pt();
713
031aac8e 714 if(fTrackType==1 || fTrackType==2)
0f76d8ae 715 track->GetImpactParametersTPC(dca2D,dcaZ); //TPConly
031aac8e 716 else
717 track->GetImpactParameters(dca2D,dcaZ); //Global
71e77a79 718
fdceab34 719 UChar_t itsMap = track->GetITSClusterMap();
fdceab34 720 for (Int_t i=0; i < 6; i++) {
721 if (itsMap & (1 << i))
722 nPointITS ++;
723 }
71e77a79 724
d889ce29 725 // fPtAll->Fill(pt);
fdceab34 726 fPtAllMC->Fill(ptMC);
71e77a79 727
fdceab34 728 if (fTrackCuts->AcceptTrack(track)) {
729
031aac8e 730 if(fTrackType==7) {
327d12da 731 if(fTrackCutsReject ) {
42881dab 732 if(fTrackCutsReject->AcceptTrack(track) ) {
733 if(track) delete track;
327d12da 734 continue;
42881dab 735 }
327d12da 736 }
737
031aac8e 738 if(esdtrack->GetConstrainedParam())
739 track->Set(esdtrack->GetConstrainedParam()->GetX(),esdtrack->GetConstrainedParam()->GetAlpha(),esdtrack->GetConstrainedParam()->GetParameter(),esdtrack->GetConstrainedParam()->GetCovariance());
031aac8e 740 }
741
742 pt = track->Pt();
743 phi = track->Phi();
744
d889ce29 745 nSigmaToVertex = fTrackCuts->GetSigmaToVertex(track);// Calculates the number of sigma to the vertex for a track.
746 chi2C = track->GetConstrainedChi2();
747 relUncertainty1Pt = TMath::Sqrt(track->GetSigma1Pt2())*pt;
748
fdceab34 749 fPtSel->Fill(pt);
36c36a0c 750 if(track->GetLabel()<0) {
751 fPtSelFakes->Fill(pt);
752 fNPointTPCFakes->Fill(track->GetTPCNcls());
753 }
fdceab34 754 fPtSelMC->Fill(ptMC);
57bee263 755 fPtAllvsPtMC->Fill(ptMC,pt);
a337a5a9 756 fPtAllminPtMCvsPtMC->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC) );
757 fPtAllminPtMCvsPtAll->Fill(pt,(1./pt-1./ptMC)/(1./ptMC) );
fdceab34 758 fPtAllminPtMCvsPtAllNPointTPC->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC),track->GetTPCNcls());
5a0bd31f 759 fPtAllminPtMCvsPtAllNPointTPCIter1->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC),track->GetTPCNclsIter1());
760 if(track->GetTPCNcls()>0.)
761 fPtAllminPtMCvsPtAllChi2TPC->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC),track->GetTPCchi2()/track->GetTPCNcls());
762 if(track->GetTPCNclsIter1()>0.)
763 fPtAllminPtMCvsPtAllChi2TPCIter1->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC),track->GetTPCchi2Iter1()/track->GetTPCNclsIter1());
fdceab34 764 fPtAllminPtMCvsPtAllDCAR->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC),dca2D);
765 fPtAllminPtMCvsPtAllDCAZ->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC),dcaZ);
766 fPtAllminPtMCvsPtAllPhi->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC),phi);
767 fPtAllminPtMCvsPtAllNPointITS->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC),nPointITS);
768 fPtAllminPtMCvsPtAllNSigmaToVertex->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC),nSigmaToVertex);
769 fPtAllminPtMCvsPtAllChi2C->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC),chi2C);
770 fPtAllminPtMCvsPtAllRel1PtUncertainty->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC),relUncertainty1Pt);
36c36a0c 771
57bee263 772 fPtAllvsPtMCvsMult->Fill(ptMC,pt,mult);
773 fPtAllminPtMCvsPtAllvsMult->Fill(ptMC,(1./pt-1./ptMC)/(1./ptMC), mult);
774
36c36a0c 775 //Check if track is reconstructed multiple times
0f76d8ae 776 /*
36c36a0c 777 int multCounter = 1;
778 for (Int_t iTrack2 = iTrack+1; iTrack2 < nTracks; iTrack2++) {
779 // AliESDtrack *track2 = GetTrackForAnalysis(iTrack2);
780 AliESDtrack *track2;
781 AliESDtrack *esdtrack2 = fESD->GetTrack(iTrack2);
782 if(!esdtrack2) continue;
783
784 if(fTrackType==1)
785 track2 = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack2->GetID());
786 else if(fTrackType==2) {
787 track2 = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack2->GetID());
788 if(!track2) {
36c36a0c 789 continue;
790 }
791 AliExternalTrackParam exParam2;
792 Bool_t relate = track2->RelateToVertexTPC(fVtx,fESD->GetMagneticField(),kVeryBig,&exParam2);
793 if( !relate ) {
794 delete track2;
795 continue;
796 }
797 track2->Set(exParam2.GetX(),exParam2.GetAlpha(),exParam2.GetParameter(),exParam2.GetCovariance());
798 }
799 else
800 track2 = esdtrack2;
801 if(!track2) {
36c36a0c 802 continue;
803 }
804
805 if (fTrackCuts->AcceptTrack(track2)) {
806 Int_t label2 = TMath::Abs(track2->GetLabel());
807 if(label==label2) {
808 fNPointTPCMultRec->Fill(track->GetTPCNcls());
809 fNPointTPCMultRec->Fill(track2->GetTPCNcls());
810 fDeltaPtMultRec->Fill(track->Pt(),track->Pt()-track2->Pt());
811 multCounter++;
812 }
813 }
814 if(fTrackType==1 || fTrackType==2) delete track2;
815 }//track2 loop
816 if(multCounter>1) fMultRec->Fill(multCounter);
0f76d8ae 817 */
36c36a0c 818
819 }//fTrackCuts selection
820
5a0bd31f 821
42881dab 822 if(fTrackType==1 || fTrackType==2 || fTrackType==7) {
5a0bd31f 823 if(track) delete track;
824 }
36c36a0c 825
fdceab34 826 }//ESD track loop
827
828 // Post output data
829 PostData(0, fHistList);
5a0bd31f 830
fdceab34 831}
b4691ee7 832//________________________________________________________________________
833Bool_t AliPWG4HighPtQAMC::PythiaInfoFromFile(const char* currFile,Float_t &fXsec,Float_t &fTrials){
834 //
835 // get the cross section and the trails either from pyxsec.root or from pysec_hists.root
836 // This is to called in Notify and should provide the path to the AOD/ESD file
837 // Copied from AliAnalysisTaskJetSpectrum2
838 //
839
840 TString file(currFile);
841 fXsec = 0;
842 fTrials = 1;
843
844 if(file.Contains("root_archive.zip#")){
845 Ssiz_t pos1 = file.Index("root_archive",12,TString::kExact);
846 Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
847 file.Replace(pos+1,20,"");
848 }
849 else {
850 // not an archive take the basename....
851 file.ReplaceAll(gSystem->BaseName(file.Data()),"");
852 }
cb76764e 853 // Printf("%s",file.Data());
b4691ee7 854
855
856 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
857 if(!fxsec){
858 // next trial fetch the histgram file
859 fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
860 if(!fxsec){
861 // not a severe condition but inciate that we have no information
862 return kFALSE;
863 }
864 else{
865 // find the tlist we want to be independtent of the name so use the Tkey
866 TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
867 if(!key){
868 fxsec->Close();
869 return kFALSE;
870 }
871 TList *list = dynamic_cast<TList*>(key->ReadObj());
872 if(!list){
873 fxsec->Close();
874 return kFALSE;
875 }
876 fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
877 fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
878 fxsec->Close();
879 }
880 } // no tree pyxsec.root
881 else {
882 TTree *xtree = (TTree*)fxsec->Get("Xsection");
883 if(!xtree){
884 fxsec->Close();
885 return kFALSE;
886 }
887 UInt_t ntrials = 0;
888 Double_t xsection = 0;
889 xtree->SetBranchAddress("xsection",&xsection);
890 xtree->SetBranchAddress("ntrials",&ntrials);
891 xtree->GetEntry(0);
892 fTrials = ntrials;
893 fXsec = xsection;
894 fxsec->Close();
895 }
896 return kTRUE;
897}
898//________________________________________________________________________
899Bool_t AliPWG4HighPtQAMC::Notify()
900{
901 //
902 // Implemented Notify() to read the cross sections
903 // and number of trials from pyxsec.root
904 // Copied from AliAnalysisTaskJetSpectrum2
905 //
906
907 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
908 Float_t xsection = 0;
909 Float_t ftrials = 1;
910
911 fAvgTrials = 1;
912 if(tree){
913 TFile *curfile = tree->GetCurrentFile();
914 if (!curfile) {
915 Error("Notify","No current file");
916 return kFALSE;
917 }
918 if(!fh1Xsec||!fh1Trials){
cb76764e 919 // Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
b4691ee7 920 return kFALSE;
921 }
922 PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
923 fh1Xsec->Fill("<#sigma>",xsection);
924 // construct a poor man average trials
925 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
926 if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
927 }
928 return kTRUE;
929}
930
931//________________________________________________________________________
932AliGenPythiaEventHeader* AliPWG4HighPtQAMC::GetPythiaEventHeader(AliMCEvent *mcEvent){
933
934 if(!mcEvent)return 0;
935 AliGenEventHeader* genHeader = mcEvent->GenEventHeader();
936 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
937 if(!pythiaGenHeader){
938 // cocktail ??
939 AliGenCocktailEventHeader* genCocktailHeader = dynamic_cast<AliGenCocktailEventHeader*>(genHeader);
940
941 if (!genCocktailHeader) {
cb76764e 942 // AliWarningGeneral(Form(" %s:%d",(char*)__FILE__,__LINE__),"Unknown header type (not Pythia or Cocktail)");
b4691ee7 943 // AliWarning(Form("%s %d: Unknown header type (not Pythia or Cocktail)",(char*)__FILE__,__LINE__));
944 return 0;
945 }
946 TList* headerList = genCocktailHeader->GetHeaders();
947 for (Int_t i=0; i<headerList->GetEntries(); i++) {
948 pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(headerList->At(i));
949 if (pythiaGenHeader)
950 break;
951 }
952 if(!pythiaGenHeader){
953 AliWarningGeneral(Form(" %s:%d",(char*)__FILE__,__LINE__),"Pythia event header not found");
954 return 0;
955 }
956 }
957 return pythiaGenHeader;
958
959}
960
fdceab34 961//________________________________________________________________________
962void AliPWG4HighPtQAMC::Terminate(Option_t *)
963{
df943115 964 // The Terminate() function is the last function to be called during
965 // a query. It always runs on the client, it can be used to present
966 // the results graphically or save the results to file.
fdceab34 967
fdceab34 968}
df943115 969
970#endif