]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/Correlations/DPhi/AliAnalysisTaskMinijet.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Adding TOF calib task for calibration of problematic channels
[u/mrichter/AliRoot.git] / PWGCF / Correlations / DPhi / AliAnalysisTaskMinijet.cxx
CommitLineData
117f99e3 1#include <TChain.h>
2#include <TList.h>
3
4#include <TTree.h>
5#include <TH1F.h>
6#include <TH2F.h>
7099c89a 7#include <THnSparse.h>
117f99e3 8#include <TProfile.h>
9#include <TCanvas.h>
10#include "TRandom.h"
11
12#include "AliVEvent.h"
13#include "AliVParticle.h"
14
15#include "AliESDEvent.h"
16#include "AliESDtrack.h"
17#include "AliMultiplicity.h"
18#include "AliESDtrackCuts.h"
19
20#include "AliAODEvent.h"
21#include "AliAODTrack.h"
b7c0438e 22#include "AliAODMCParticle.h"
7099c89a 23#include "AliAODMCHeader.h"
117f99e3 24
25#include "AliStack.h"
26#include "AliMCEvent.h"
27#include "AliMCParticle.h"
28#include "AliGenEventHeader.h"
fa8e4941 29#include "AliCentrality.h"
117f99e3 30
31#include "AliAnalysisManager.h"
7099c89a 32#include "AliInputEventHandler.h"
33
34#include <vector>
35#include <algorithm>
36#include <functional>
37using namespace std;
117f99e3 38
39#include "AliAnalysisTaskMinijet.h"
40
f2565a6b 41// Analysis task for two-particle correlations using all particles over pt threshold
42// pt_trig threshold for trigger particle (event axis) and pt_assoc for possible associated particles.
43// Extract mini-jet yield and fragmentation properties via Delta-Phi histograms of these correlations
44// post processing of analysis output via macro plot3and2Gaus.C
45// Can use ESD or AOD, reconstructed and Monte Carlo data as input
fa8e4941 46// Author: Eva Sicking, modifications by Emilia Leogrande
f2565a6b 47
117f99e3 48
49ClassImp(AliAnalysisTaskMinijet)
50
51//________________________________________________________________________
fa8e4941 52AliAnalysisTaskMinijet::AliAnalysisTaskMinijet(const char *name)
b9ad6f04 53 : AliAnalysisTaskSE(name),
fa8e4941 54 fUseMC(kFALSE),
55 fMcOnly(kFALSE),
56 fBSign(0),
57 fAnalysePrimOnly(kFALSE),// not used
58 fPtMin(0.2),
59 fPtMax(50.0),
60 fCuts(0),
61 fTriggerPtCut(0.7),
62 fAssociatePtCut(0.4),
63 fMode(0),
64 fTriggerType(1),
65 fFilterBit(128),
66 fVertexZCut(10.),
67 fEtaCut(0.9),
68 fEtaCutSeed(0.9),
69 fSelectParticles(1),
70 fSelectParticlesAssoc(1),
71 fCheckSDD(true),
72 fSelOption(1),
73 fCorrStrangeness(true),
74 fThreeParticleCorr(false),
75 fRejectChunks(false),
76 fNTPC(10),
77 fESDEvent(0),
78 fAODEvent(0),
79 fNMcPrimAccept(0),
80 fNRecAccept(0),
81 fNRecAcceptStrangeCorr(0),
82 fNMcPrimAcceptTracklet(0),
83 fNRecAcceptTracklet(0),
84 fVzEvent(0),
85 fMeanPtRec(0),
86 fLeadingPtRec(0),
87 fHists(0),
88 fStep(0),
89 fEventStat(0),
90 fHistPt(0),
91 fHistPtMC(0),
92 fNContrNtracklets(0),
93 fNContrNtracks(0),
94 fCorruptedChunks(0),
95 fCorruptedChunksAfter(0),
96 fNmcNch(0),
97 fPNmcNch(0),
98 fNmcNchVtx(0),
99 fNmcNchVtxStrangeCorr(0),
100 fPNmcNchVtx(0),
101 fNmcNchTracklet(0),
102 fPNmcNchTracklet(0),
103 fNmcNchVtxTracklet(0),
104 fPNmcNchVtxTracklet(0),
105 fChargedPi0(0),
106 fVertexCheck(0),
107 fPropagateDca(0),
3eb806f8 108 fCentralityMethod("")
117f99e3 109{
fa8e4941 110
111 //Constructor
112
113 for(Int_t i = 0;i< 8;i++){
114 fMapSingleTrig[i] = 0;
115 fMapPair[i] = 0;
116 fMapEvent[i] = 0;
117 fMapAll[i] = 0;
118 fMapThree[i] = 0;
119
120 fVertexZ[i] = 0;
121
122 fNcharge[i] = 0;
123 fPt[i] = 0;
124 fEta[i] = 0;
125 fPhi[i] = 0;
126 fDcaXY[i] = 0;
127 fDcaZ[i] = 0;
128
129 fPtSeed[i] = 0;
130 fEtaSeed[i] = 0;
131 fPhiSeed[i] = 0;
132
133 fPtOthers[i] = 0;
134 fEtaOthers[i] = 0;
135 fPhiOthers[i] = 0;
136 fPtEtaOthers[i] = 0;
137
138 fPhiEta[i] = 0;
139
140 fDPhiDEtaEventAxis[i] = 0;
141 fDPhiDEtaEventAxisSeeds[i]= 0;
142 fTriggerNch[i] = 0;
143
144 fTriggerNchSeeds[i] = 0;
145 fTriggerTracklet[i] = 0;
146
147 fNch07Nch[i] = 0;
148 fPNch07Nch[i] = 0;
149
150 fNch07Tracklet[i] = 0;
151 fNchTracklet[i] = 0;
152 fPNch07Tracklet[i] = 0;
153 fDPhiEventAxis[i] = 0;
154 fDPhi1DPhi2[i] = 0;
155 }
156 DefineOutput(1, TList::Class());
117f99e3 157}
158
159//________________________________________________________________________
160AliAnalysisTaskMinijet::~AliAnalysisTaskMinijet()
161{
fa8e4941 162 // Destructor
163
164 if (fHists && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fHists;
117f99e3 165}
166
167//________________________________________________________________________
168void AliAnalysisTaskMinijet::UserCreateOutputObjects()
169{
fa8e4941 170 // Create histograms
171 // Called once
172 if(fDebug) Printf("In User Create Output Objects.");
173
7d6dc277 174 Int_t nbinsCentr = 0;
175 Float_t minbinCentr=0, maxbinCentr=0;
176
177 if (fCentralityMethod.Length() > 0)
178 {
0a24e635 179 nbinsCentr = 105;
7d6dc277 180 minbinCentr=0;
bf937b73 181 maxbinCentr=105;
7d6dc277 182 }
183 else
184 {
185 nbinsCentr = 101;
186 minbinCentr=-0.5;
187 maxbinCentr=100.5;
188 }
fa8e4941 189
190 fStep = new TH1F("fStep", "fStep", 10, -0.5, 9.5);
191 fEventStat = new TH1F("fEventStat", "fEventStat", 10, -0.5, 9.5);
192 fHistPt = new TH1F("fHistPt", "P_{T} distribution REC", 150, 0.1, 3.1);
193 fHistPt->GetXaxis()->SetTitle("P_{T} (GeV/c)");
194 fHistPt->GetYaxis()->SetTitle("dN/dP_{T} (c/GeV)");
195 fHistPt->SetMarkerStyle(kFullCircle);
196 fNContrNtracklets = new TH2F ("fNContrNtracklets", ";N_{tracklets};N_{vtx contrib}", nbinsCentr,minbinCentr,maxbinCentr,nbinsCentr,minbinCentr,maxbinCentr);
197 fNContrNtracks = new TH2F ("fNContrNtracks", ";N_{tracks};N_{vtx contrib}", nbinsCentr,minbinCentr,maxbinCentr,nbinsCentr,minbinCentr,maxbinCentr);
198 fCorruptedChunks = new TH2F ("fCorruptedChunks",
199 ";N_{tracks,TPC};N_{tracks,ITS-TPC}", nbinsCentr,minbinCentr,maxbinCentr,nbinsCentr,minbinCentr,maxbinCentr);
200 fCorruptedChunksAfter = new TH2F ("fCorruptedChunksAfter",
201 ";N_{tracks,TPC};N_{tracks,ITS-TPC}", nbinsCentr,minbinCentr,maxbinCentr,nbinsCentr,minbinCentr,maxbinCentr);
7099c89a 202
117f99e3 203
fa8e4941 204
205
206 if (fUseMC) {
207 fHistPtMC = new TH1F("fHistPtMC", "P_{T} distribution MC", 150, 0.1, 3.1);
208 fHistPtMC->GetXaxis()->SetTitle("P_{T} (GeV/c)");
209 fHistPtMC->GetYaxis()->SetTitle("dN/dP_{T} (c/GeV)");
210 fHistPtMC->SetMarkerStyle(kFullCircle);
211
212 fNmcNch = new TH2F("fNmcNch", "fNmcNch", nbinsCentr,minbinCentr,maxbinCentr,nbinsCentr,minbinCentr,maxbinCentr);
213 fPNmcNch = new TProfile("pNmcNch", "pNmcNch", nbinsCentr,minbinCentr,maxbinCentr);
214 fNmcNchVtx = new TH2F("fNmcNchVtx", "fNmcNchVtx", nbinsCentr,minbinCentr,maxbinCentr,nbinsCentr,minbinCentr,maxbinCentr);
215 fNmcNchVtxStrangeCorr = new TH2F("fNmcNchVtxStrangeCorr", "fNmcNchVtxStrangeCorr", nbinsCentr,minbinCentr,maxbinCentr,nbinsCentr,minbinCentr,maxbinCentr);
216 fPNmcNchVtx = new TProfile("pNmcNchVtx", "pNmcNchVtx", nbinsCentr,minbinCentr,maxbinCentr);
217
218 fNmcNchTracklet = new TH2F("fNmcNchTracklet", "fNmcNchTracklet", nbinsCentr,minbinCentr,maxbinCentr,nbinsCentr,minbinCentr,maxbinCentr);
219 fPNmcNchTracklet = new TProfile("pNmcNchTracklet", "pNmcNchTracklet", nbinsCentr,minbinCentr,maxbinCentr);
220 fNmcNchVtxTracklet = new TH2F("fNmcNchVtxTracklet", "fNmcNchVtxTracklet", nbinsCentr,minbinCentr,maxbinCentr,nbinsCentr,minbinCentr,maxbinCentr);
221 fPNmcNchVtxTracklet = new TProfile("pNmcNchVtxTracklet", "pNmcNchVtxTracklet", nbinsCentr,minbinCentr,maxbinCentr);
222
223
224
225 }
226
227 fChargedPi0 = new TH2F("fChargedPi0", "fChargedPi0", 200, -0.5, 199.5, 200, -0.5, 199.5);
228 fVertexCheck = new TH1F("fVertexCheck", "fVertexCheck", 100, -0.5, 99.5);
229 fPropagateDca = new TH1F("fPropagateDca", "fPropagateDca", 2, -0.5, 1.5);
230
231 //----------------------
232 //bins for pt in THnSpare
233 Double_t ptMin = 0.0, ptMax = 100.;
234 Int_t nPtBins = 37;
235 Double_t binsPt[] = {0.0, 0.1, 0.2, 0.3, 0.35, 0.4, 0.45, 0.5, 0.6, 0.7, 0.8,
236 0.9, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 6.0, 7.0, 8.0, 9.0,
237 10.0, 12.0, 14.0, 16.0, 18.0, 20.0, 25.0, 30.0, 35.0, 40.0, 45.0, 50.0, 100.0};
238
239 //Int_t nPtBinsAll = 100;
240 //Double_t ptMinAll = 1.e-2, ptMaxAll = 100.;
241 //Double_t *binsPtAll = 0;
242 //binsPtAll = (Double_t*)CreateLogAxis(nPtBinsAll,ptMinAll,ptMaxAll);
243
244
245 // Double_t ptMin2 = 0.0, ptMax2 = 100.;
246 // Int_t nPtBins2 = 9;
247 // Double_t binsPt2[] = {0.1, 0.4, 0.7, 1.0, 2.0, 3.0, 4.0, 5.0, 10.0, 100.0};
248
249 // Int_t nPtBins2 = 10;
250 // Double_t ptMin2 = 0.4, ptMax2 = 100.;
251 // Double_t *binsPt2 = 0;
252 // binsPt2 = CreateLogAxis(nPtBins2,ptMin2,ptMax2);
253
254 //3 dim matrix
255 Int_t binsEffHisto[3] = {Int_t(fEtaCut*20), nPtBins, nbinsCentr };
256 Double_t minEffHisto[3] = {-fEtaCut, ptMin, minbinCentr };
257 Double_t maxEffHisto[3] = {fEtaCut, ptMax, maxbinCentr};
258
259 //5 dim matrix
dd37538a 260 Int_t binsEffHisto5[6] = { nPtBins, nPtBins, 1, 90, nbinsCentr , 2 };
fa8e4941 261 Double_t minEffHisto5[6] = { ptMin, ptMin, -2*fEtaCut, -0.5*TMath::Pi(), minbinCentr , -0.5 };
262 Double_t maxEffHisto5[6] = { ptMax, ptMax, 2*fEtaCut, 1.5*TMath::Pi(), maxbinCentr , 1.5 };
263
264
265 //4 dim matrix
266 Int_t binsEvent[4] = { nbinsCentr, 20, 50, nPtBins };
267 Double_t minEvent[4] = { minbinCentr, -10, 0, ptMin };
268 Double_t maxEvent[4] = { maxbinCentr, 10, 10, ptMax };
269
270 //3 dim matrix
271 Int_t binsAll[3] = {Int_t(fEtaCut*20), nPtBins, nbinsCentr };
272 Double_t minAll[3] = {-fEtaCut, ptMin, minbinCentr };
273 Double_t maxAll[3] = {fEtaCut, ptMax, maxbinCentr };
274
275 //3 dim matrix
276 Int_t binsThree[3] = { 90, 90, nbinsCentr};
277 Double_t minThree[3] = {-0.5*TMath::Pi(), -0.5*TMath::Pi(), minbinCentr};
278 Double_t maxThree[3] = { 1.5*TMath::Pi(), 1.5*TMath::Pi(), maxbinCentr};
279
7099c89a 280
fa8e4941 281
282
283
284 //--------------------
285 TString dataType[2] ={"ESD", "AOD"};
286 TString labels[8]={Form("%sAllAllMcNmc",dataType[fMode].Data()),
287 Form("%sTrigAllMcNmc",dataType[fMode].Data()),
288 Form("%sTrigVtxMcNmc",dataType[fMode].Data()),
289 Form("%sTrigVtxMcNrec",dataType[fMode].Data()),
290 Form("%sTrigVtxRecMcPropNrec",dataType[fMode].Data()),
291 Form("%sTrigVtxRecNrec",dataType[fMode].Data()),
292 Form("%sTrigVtxRecMcPropNrecStrangeCorr",dataType[fMode].Data()),
293 Form("%sTrigVtxRecNrecStrangeCorr",dataType[fMode].Data())};
294
295
296 for(Int_t i=0;i<8;i++){
297
298 fMapSingleTrig[i] = new THnSparseD(Form("fMapSingleTrig%s", labels[i].Data()),"eta:pt:Nrec",
299 3,binsEffHisto,minEffHisto,maxEffHisto);
300 fMapSingleTrig[i]->SetBinEdges(1,binsPt);
301 fMapSingleTrig[i]->GetAxis(0)->SetTitle("#eta");
302 fMapSingleTrig[i]->GetAxis(1)->SetTitle("p_{T} (GeV/c)");
303 fMapSingleTrig[i]->GetAxis(2)->SetTitle("N_{rec}");
304 fMapSingleTrig[i]->Sumw2();
305
306 fMapPair[i] = new THnSparseD(Form("fMapPair%s", labels[i].Data()),"pt_trig:pt_assoc:DeltaEta:DeltaPhi:Nrec:likesign",
307 6,binsEffHisto5,minEffHisto5,maxEffHisto5);
308 fMapPair[i]->SetBinEdges(0,binsPt);
309 fMapPair[i]->SetBinEdges(1,binsPt);
310 fMapPair[i]->GetAxis(0)->SetTitle("p_{T, trig} (GeV/c)");
311 fMapPair[i]->GetAxis(1)->SetTitle("p_{T, assoc} (GeV/c)");
312 fMapPair[i]->GetAxis(2)->SetTitle("#Delta #eta");
313 fMapPair[i]->GetAxis(3)->SetTitle("#Delta #phi");
314 fMapPair[i]->GetAxis(4)->SetTitle("N_{rec}");
315 fMapPair[i]->GetAxis(5)->SetTitle("Like-sign or Unlike-sign");
316 fMapPair[i]->Sumw2();
317
318
319 fMapEvent[i] = new THnSparseD(Form("fMapEvent%s", labels[i].Data()),"Nrec:vertexZ:meanPt:leadingPt",
320 4,binsEvent,minEvent,maxEvent);
321 fMapEvent[i]->GetAxis(0)->SetTitle("N_{rec}");
322 fMapEvent[i]->GetAxis(1)->SetTitle("z_{vertex} (cm)");
323 fMapEvent[i]->GetAxis(2)->SetTitle("meanPt");
324 fMapEvent[i]->SetBinEdges(3,binsPt);
325 fMapEvent[i]->GetAxis(3)->SetTitle("leadingPt");
326 fMapEvent[i]->Sumw2();
327
328 fMapAll[i] = new THnSparseD(Form("fMapAll%s", labels[i].Data()),"eta:pt:Nrec",
329 3,binsAll,minAll,maxAll);
330 fMapAll[i]->SetBinEdges(1,binsPt);
331 fMapAll[i]->GetAxis(0)->SetTitle("#eta");
332 fMapAll[i]->GetAxis(1)->SetTitle("p_{T} (GeV/c)");
333 fMapAll[i]->GetAxis(2)->SetTitle("N_{rec}");
334 fMapAll[i]->Sumw2();
335
336
337 fMapThree[i] = new THnSparseD(Form("fMapThree%s", labels[i].Data()),"dphi1:dphi2:Nrec",
338 3,binsThree,minThree,maxThree);
339 fMapThree[i]->GetAxis(0)->SetTitle("#Delta#varphi_{1}");
340 fMapThree[i]->GetAxis(1)->SetTitle("#Delta#varphi_{2}");
341 fMapThree[i]->GetAxis(2)->SetTitle("N_{rec}");
342 fMapThree[i]->Sumw2();
343
344
345 fVertexZ[i] = new TH1F(Form("fVertexZ%s",labels[i].Data()),
346 Form("fVertexZ%s",labels[i].Data()) ,
347 220, -11., 11.);
348 fPt[i] = new TH1F(Form("fPt%s",labels[i].Data()),
349 Form("fPt%s",labels[i].Data()) ,
350 200, 0., 50);
351 fNcharge[i] = new TH1F(Form("fNcharge%s",labels[i].Data()),
352 Form("fNcharge%s",labels[i].Data()) ,
353 250, -0.5, 249.5);
354 fEta[i] = new TH1F (Form("fEta%s",labels[i].Data()),
355 Form("fEta%s",labels[i].Data()) ,
356 100, -2., 2);
357 fPhi[i] = new TH1F(Form("fPhi%s",labels[i].Data()),
358 Form("fPhi%s",labels[i].Data()) ,
359 360, 0.,2*TMath::Pi());
360 fDcaXY[i] = new TH1F(Form("fDcaXY%s",labels[i].Data()),
361 Form("fDcaXY%s",labels[i].Data()) ,
362 200, -0.3,0.3);
363 fDcaZ[i] = new TH1F(Form("fDcaZ%s",labels[i].Data()),
364 Form("fDcaZ%s",labels[i].Data()) ,
365 200, -2.2,2.2);
366 fPtSeed[i] = new TH1F(Form("fPSeedt%s",labels[i].Data()),
367 Form("fPtSeed%s",labels[i].Data()) ,
368 500, 0., 50);
369 fEtaSeed[i] = new TH1F (Form("fEtaSeed%s",labels[i].Data()),
370 Form("fEtaSeed%s",labels[i].Data()) ,
371 100, -2., 2);
372 fPhiSeed[i] = new TH1F(Form("fPhiSeed%s",labels[i].Data()),
373 Form("fPhiSeed%s",labels[i].Data()) ,
374 360, 0.,2*TMath::Pi());
375
376 fPtOthers[i] = new TH1F(Form("fPOtherst%s",labels[i].Data()),
377 Form("fPtOthers%s",labels[i].Data()) ,
378 500, 0., 50);
379 fEtaOthers[i] = new TH1F (Form("fEtaOthers%s",labels[i].Data()),
380 Form("fEtaOthers%s",labels[i].Data()) ,
381 100, -2., 2);
382 fPhiOthers[i] = new TH1F(Form("fPhiOthers%s",labels[i].Data()),
383 Form("fPhiOthers%s",labels[i].Data()) ,
384 360, 0.,2*TMath::Pi());
385 fPtEtaOthers[i] = new TH2F(Form("fPtEtaOthers%s",labels[i].Data()),
386 Form("fPtEtaOthers%s",labels[i].Data()) ,
387 500, 0., 50, 100, -1., 1);
388 fPhiEta[i] = new TH2F(Form("fPhiEta%s",labels[i].Data()),
389 Form("fPhiEta%s",labels[i].Data()) ,
390 180, 0., 2*TMath::Pi(), 100, -1.,1.);
391 fDPhiDEtaEventAxis[i] = new TH2F(Form("fDPhiDEtaEventAxis%s",labels[i].Data()),
392 Form("fDPhiDEtaEventAxis%s",labels[i].Data()) ,
393 180, -0.5* TMath::Pi(), 1.5*TMath::Pi(), 9, -1.8,1.8);
394 fTriggerNch[i] = new TH1F(Form("fTriggerNch%s",labels[i].Data()),
395 Form("fTriggerNch%s",labels[i].Data()) ,
396 250, -0.5, 249.5);
397 fTriggerNchSeeds[i] = new TH2F(Form("fTriggerNchSeeds%s",labels[i].Data()),
398 Form("fTriggerNchSeeds%s",labels[i].Data()) ,
399 250, -0.5, 249.5, 100, -0.5, 99.5);
400 fTriggerTracklet[i] = new TH1F(Form("fTriggerTracklet%s",labels[i].Data()),
401 Form("fTriggerTracklet%s",labels[i].Data()) ,
402 250, -0.5, 249.5);
403 fNch07Nch[i] = new TH2F(Form("fNch07Nch%s",labels[i].Data()),
404 Form("fNch07Nch%s",labels[i].Data()) ,
405 250, -2.5, 247.5,250, -2.5, 247.5);
406 fPNch07Nch[i] = new TProfile(Form("pNch07Nch%s",labels[i].Data()),
407 Form("pNch07Nch%s",labels[i].Data()) ,
408 250, -2.5, 247.5);
409 fNch07Tracklet[i] = new TH2F(Form("fNch07Tracklet%s",labels[i].Data()),
410 Form("fNch07Tracklet%s",labels[i].Data()) ,
411 250, -2.5, 247.5,250, -2.5, 247.5);
412 fNchTracklet[i] = new TH2F(Form("fNchTracklet%s",labels[i].Data()),
413 Form("fNchTracklet%s",labels[i].Data()) ,
414 250, -2.5, 247.5,250, -2.5, 247.5);
415 fPNch07Tracklet[i] = new TProfile(Form("pNch07Tracklet%s",labels[i].Data()),
416 Form("pNch07Tracklet%s",labels[i].Data()) ,
417 250, -2.5, 247.5);
418 fDPhiEventAxis[i] = new TH1F(Form("fDPhiEventAxis%s",
419 labels[i].Data()),
420 Form("fDPhiEventAxis%s",
421 labels[i].Data()) ,
422 180, -0.5*TMath::Pi(), 1.5*TMath::Pi());
423 fDPhi1DPhi2[i] = new TH2F(Form("fDPhi1DPhi2%s",labels[i].Data()),
424 Form("fDPhi1DPhi2%s",labels[i].Data()),
425 180, -0.5* TMath::Pi(), 1.5*TMath::Pi(),
426 180, -0.5* TMath::Pi(), 1.5*TMath::Pi());
427 }
428
429 fHists = new TList();
430 fHists->SetOwner();
431
432 fHists->Add(fStep);
433 fHists->Add(fEventStat);
434 fHists->Add(fHistPt);
435 fHists->Add(fNContrNtracklets);
436 fHists->Add(fNContrNtracks);
437 fHists->Add(fCorruptedChunks);
438 fHists->Add(fCorruptedChunksAfter);
439
440 if(fUseMC){
441 fHists->Add(fHistPtMC);
442
443 fHists->Add(fNmcNch);
444 fHists->Add(fPNmcNch);
445 fHists->Add(fNmcNchVtx);
446 fHists->Add(fNmcNchVtxStrangeCorr);
447 fHists->Add(fPNmcNchVtx);
448
449 fHists->Add(fNmcNchTracklet);
450 fHists->Add(fPNmcNchTracklet);
451 fHists->Add(fNmcNchVtxTracklet);
452 fHists->Add(fPNmcNchVtxTracklet);
453 }
454 fHists->Add(fChargedPi0);
455 fHists->Add(fVertexCheck);
456 fHists->Add(fPropagateDca);
457
458 for(Int_t i=0;i<8;i++){
459 fHists->Add(fMapSingleTrig[i]);
460 fHists->Add(fMapPair[i]);
461 fHists->Add(fMapEvent[i]);
462 fHists->Add(fMapAll[i]);
463 fHists->Add(fMapThree[i]);
464 fHists->Add(fVertexZ[i]);
465 fHists->Add(fPt[i]);
466 fHists->Add(fNcharge[i]);
467 fHists->Add(fEta[i]);
468 fHists->Add(fPhi[i]);
469 fHists->Add(fDcaXY[i]);
470 fHists->Add(fDcaZ[i]);
471 fHists->Add(fPtSeed[i]);
472 fHists->Add(fEtaSeed[i]);
473 fHists->Add(fPhiSeed[i]);
474 fHists->Add(fPtOthers[i]);
475 fHists->Add(fEtaOthers[i]);
476 fHists->Add(fPhiOthers[i]);
477 fHists->Add(fPtEtaOthers[i]);
478 fHists->Add(fPhiEta[i]);
479 fHists->Add(fDPhiDEtaEventAxis[i]);
480 fHists->Add(fTriggerNch[i]);
481 fHists->Add(fTriggerNchSeeds[i]);
482 fHists->Add(fTriggerTracklet[i]);
483 fHists->Add(fNch07Nch[i]);
484 fHists->Add(fPNch07Nch[i]);
485 fHists->Add(fNch07Tracklet[i]);
486 fHists->Add(fNchTracklet[i]);
487 fHists->Add(fPNch07Tracklet[i]);
488 fHists->Add(fDPhiEventAxis[i]);
489 fHists->Add(fDPhi1DPhi2[i]);
490 }
491
492 PostData(1, fHists);
493
117f99e3 494}
495
496//________________________________________________________________________
497void AliAnalysisTaskMinijet::UserExec(Option_t *)
498{
fa8e4941 499 // Main function, called for each event
500 // Kinematics-only, ESD and AOD can be processed.
501 // Data is read (ReadEventESD, ReadEventAOD...) and then analysed (Analyse).
502 // - in case of MC with full detector simulation, all correction steps(0-5) can be processed
503 // - for Data, only step 5 is performed
504 // - for kinematics-only, only step 0 is processed
505
506 // - trigger - - vertex - - tracks - - multiplicity -
507 // step 7 = Triggered events, reconstructed accepted vertex, reconstructed tracks + strangness corr, reconstructed multiplicity+strangeCorr
508 // step 6 = Triggered events, reconstructed accepted vertex, reconstructed tracks with MC prop + stangess corr, reconstructed multiplicity+strangeCorr
509 // step 5 = Triggered events, reconstructed accepted vertex, reconstructed tracks, reconstructed multiplicity
510 // step 4 = Triggered events, reconstructed accepted vertex, reconstructed tracks with MC properties, reconstructed multiplicity
511 // step 3 = Triggered events, reconstructed accepted vertex, mc primary particles, reconstructed multiplicity
512 // step 2 = Triggered events, reconstructed accepted vertex, mc primary particles, true multiplicity
513 // step 1 = Triggered events, all mc primary particles, true multiplicity
514 // step 0 = All events, all mc primary particles, true multiplicity
515
516
517 if(fDebug) Printf("UserExec: Event starts");
518
519 AliVEvent *event = InputEvent();
520 if (!event) {
521 Error("UserExec", "Could not retrieve event");
522 return;
ac305649 523 }
fa8e4941 524 fBSign= event->GetMagneticField();
525
526 //get events, either ESD or AOD
527 fESDEvent = dynamic_cast<AliESDEvent*> (InputEvent());
528 fAODEvent = dynamic_cast<AliAODEvent*> (InputEvent());
529
530 vector<Float_t> pt;
531 vector<Float_t> eta;
532 vector<Float_t> phi;
533 vector<Short_t> charge;
534 vector<Float_t> strangenessWeight;
535 vector<Float_t> px;
536 vector<Float_t> py;
537 vector<Float_t> pz;
538 vector<Float_t> theta;
539
540
541 //number of accepted tracks and tracklets
11722f75 542 Double_t ntracks = 0; //return value for reading functions for ESD and AOD
fa8e4941 543 //Int_t ntracksRemove = 0; //return value for reading functions for ESD and AOD
11722f75 544 vector<Double_t> nTracksTracklets; // [0]=nAccepted,1=ntracklets,2=nall(also neutral in case of mc,
fa8e4941 545 //for real nall=ncharged)
546
547 if(!fAODEvent && !fESDEvent)return;
548
549 //Centrality
550 Double_t centrality = 0;
551 if (fCentralityMethod.Length() > 0)
552 {
553 AliCentrality *centralityObj = 0;
554 if (fAODEvent)
555 centralityObj = fAODEvent->GetHeader()->GetCentralityP();
556 else if (fESDEvent)
557 centralityObj = fESDEvent->GetCentrality();
558 if (centralityObj)
559 {
560 centrality = centralityObj->GetCentralityPercentile(fCentralityMethod);
561 AliInfo(Form("Centrality is %f", centrality));
562 }
563 else
564 {
565 Printf("WARNING: Centrality object is 0");
566 centrality = -1;
567 }
ac305649 568 }
fa8e4941 569
570 // SDD check for LHC11a
571 if (fCheckSDD) {
572
573 //ESD
574 if(fESDEvent){
575 if(fSelOption==0){
576 const AliMultiplicity *mul = fESDEvent->GetMultiplicity();
577 Int_t nClu3 = mul->GetNumberOfITSClusters(2);
578 Int_t nClu4 = mul->GetNumberOfITSClusters(3);
579 if(nClu3==0 && nClu4==0) return;
580 }
581 else if (fSelOption==1){
582 TString trcl = fESDEvent->GetFiredTriggerClasses().Data();
583 if (!(trcl.Contains("CINT1-B-NOPF-ALLNOTRD"))) return;
584 }
585 }
586
587 //AOD
588 if(fAODEvent){
589 if(fSelOption==0){
590 Bool_t useEvent = false;
591 Int_t nTracks = fAODEvent->GetNTracks();
592 for(Int_t itrack=0; itrack<nTracks; itrack++) {
593 AliAODTrack * track = fAODEvent->GetTrack(itrack);
594 if(TESTBIT(track->GetITSClusterMap(),2) || TESTBIT(track->GetITSClusterMap(),3) ){
595 useEvent=true;
596 break;
597 }
598 }
599 if (!useEvent) return;
600 }
601 else if(fSelOption==1){
602 TString trcl = fAODEvent->GetFiredTriggerClasses().Data();
603 if (!(trcl.Contains("CINT1-B-NOPF-ALLNOTRD"))) return;
604 }
605 }
8855d81c 606 }
fa8e4941 607
608 //reset values
609 fNMcPrimAccept=0;// number of accepted primaries
610 fNRecAccept=0; // number of accepted tracks
611 fNRecAcceptStrangeCorr=0; // number of accepted tracks + strangeness correction
612 fNMcPrimAcceptTracklet=0;// number of accepted primaries (no pt cut)
613 fNRecAcceptTracklet=0; // number of accepted tracklets
614
615 // instead of task->SelectCollisionCandidate(mask) in AddTask macro
616 Bool_t isSelected = (((AliInputEventHandler*)
617 (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))
618 ->IsEventSelected() & fTriggerType);
619
620
621 if(fDebug){
622 Printf("IsSelected = %d", isSelected);
623 Printf("CheckEvent(true)= %d", CheckEvent(true));
624 Printf("CheckEvent(false)= %d", CheckEvent(false));
117f99e3 625 }
fa8e4941 626
627 fEventStat->Fill(0);//all events
628
629 //check trigger
630 if(isSelected){ // has offline trigger
631
632 fEventStat->Fill(1);//triggered event
633
634 if(CheckEvent(true)){//step 5 = TrigVtxRecNrec, step 4 = TrigVtxRecMcPropNrec ,step 3 = TrigVtxMcNrec
635
636 fEventStat->Fill(2);//triggered event with vertex
637
638 if(!fMcOnly){
639 //step 5 = TrigVtxRecNrec
640
641 // read tracks
642 if(fESDEvent) ntracks = ReadEventESD(pt, eta, phi, charge,strangenessWeight, nTracksTracklets, 5);
643 else if(fAODEvent)ntracks = ReadEventAOD(pt, eta, phi, charge,strangenessWeight, nTracksTracklets, 5);
28cc9dd5 644 else AliInfo("Fatal Error");
fa8e4941 645
646 if (fCentralityMethod.Length() > 0)
11722f75 647 ntracks = centrality;
fa8e4941 648
649 // analyse
650 if(pt.size()){ //(internally ntracks=fNRecAccept)
651 fEventStat->Fill(3);//triggered event with vertex and one reconstructed track in acceptance
652 Analyse(pt, eta, phi, charge, strangenessWeight, ntracks, nTracksTracklets[1], nTracksTracklets[2], 5);//step 5 = TrigVtxRecNrec
653 }
654
655 if(fCorrStrangeness){
656 //step 7 = TrigVtxRecNrecStrangeCorr
657
658 // read tracks
659 if(fESDEvent) ntracks = ReadEventESD(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 7);//stagness version not yet implemented
660 else if(fAODEvent)ntracks = ReadEventAOD(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 7);
28cc9dd5 661 else AliInfo("Fatal Error");
fa8e4941 662
3eb806f8 663 if (fCentralityMethod.Length() > 0)
11722f75 664 ntracks = centrality;
3eb806f8 665
fa8e4941 666 // analyse
667 if(pt.size()){ //(internally ntracks=fNRecAccept)
11722f75 668 Analyse(pt, eta, phi, charge, strangenessWeight, fNRecAcceptStrangeCorr, nTracksTracklets[1], nTracksTracklets[2], 7);//step 7 = TrigVtxRecNrecStrangeCorr
fa8e4941 669 }
670 }
671
672 if(fUseMC){
673 // step 4 = TrigVtxRecMcPropNrec
674
675 // read tracks
676 if(fESDEvent) ntracks = ReadEventESDRecMcProp(pt, eta, phi, charge,strangenessWeight, nTracksTracklets, 4);
677 else if(fAODEvent) ntracks = ReadEventAODRecMcProp(pt, eta, phi, charge,strangenessWeight, nTracksTracklets, 4);
28cc9dd5 678 else AliInfo("Fatal Error");
fa8e4941 679
3eb806f8 680 if (fCentralityMethod.Length() > 0)
a2ead254 681 ntracks = centrality;
3eb806f8 682
fa8e4941 683 //analyse
684 if(pt.size()){//(internally ntracks=fNRecAccept)
685 Analyse(pt, eta, phi, charge,strangenessWeight, ntracks, nTracksTracklets[1], nTracksTracklets[2], 4); //step 4 = TrigVtxRecMcPropNrec
686 }
687
688
689 if(fCorrStrangeness){
690 // step 6 = TrigVtxRecMcPropNrecStrangeCorr
691
692 // read tracks
693 if(fESDEvent) ntracks = ReadEventESDRecMcProp(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 6);//stagness version not yet implemented
694 else if(fAODEvent) ntracks = ReadEventAODRecMcProp(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 6);
28cc9dd5 695 else AliInfo("Fatal Error");
fa8e4941 696
3eb806f8 697 if (fCentralityMethod.Length() > 0)
11722f75 698 ntracks = centrality;
3eb806f8 699
fa8e4941 700 //analyse
701 if(pt.size()){//(internally ntracks=fNRecAccept)
11722f75 702 Analyse(pt, eta, phi, charge, strangenessWeight, fNRecAcceptStrangeCorr, nTracksTracklets[1], nTracksTracklets[2], 6); //step 6 = TrigVtxRecMcPropNrecStrangeCorr
fa8e4941 703 }
704 }
705 // step 3 = TrigVtxMcNrec
706
707 // read tracks
708 if(fESDEvent) ntracks = ReadEventESDMC(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 3);
709 else if(fAODEvent) ntracks = ReadEventAODMC(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 3);
28cc9dd5 710 else AliInfo("Fatal Error");
fa8e4941 711
3eb806f8 712 if (fCentralityMethod.Length() > 0){
11722f75 713 fNRecAccept = centrality;
714 fNMcPrimAccept = centrality;
3eb806f8 715 }
716
fa8e4941 717 // analyse
718 if(pt.size()){
719 Analyse(pt, eta, phi, charge, strangenessWeight, fNRecAccept, nTracksTracklets[1],nTracksTracklets[2], 3); //step 3 = TrigVtxMcNrec
720 Analyse(pt, eta, phi, charge, strangenessWeight, fNMcPrimAccept, nTracksTracklets[1],nTracksTracklets[2], 2); //step 2 = TrigVtxMcNmc
721 }
722
723 }
724 }
725
726 }//check event (true)
727
728
729 if(fUseMC && !fMcOnly){
730 //reset values
731 fNMcPrimAccept=0;// number of accepted primaries
732 fNRecAccept=0; // number of accepted tracks
733 fNMcPrimAcceptTracklet=0;// number of accepted primaries (no pt cut)
734 fNRecAcceptTracklet=0; // number of accepted tracklets
735
736 if(CheckEvent(false)){// all events, with and without reconstucted vertex
737 if(fESDEvent) ntracks = ReadEventESDMC(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 1);//read tracks
738 else if(fAODEvent) ntracks = ReadEventAODMC(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 1);//read tracks
28cc9dd5 739 else AliInfo("Fatal Error");
fa8e4941 740
3eb806f8 741 if (fCentralityMethod.Length() > 0)
11722f75 742 fNMcPrimAccept = centrality;
3eb806f8 743
744
fa8e4941 745 // analyse
746 if(pt.size()){
747 Analyse(pt, eta, phi, charge, strangenessWeight, fNMcPrimAccept, nTracksTracklets[1],nTracksTracklets[2], 1); // step 1 = TrigAllMcNmc
748
749 Analyse(pt, eta, phi, charge, strangenessWeight, fNMcPrimAccept, nTracksTracklets[1],nTracksTracklets[2], 0); //first part of step 0 // step 0 = AllAllMcNmc
750 }
751
752
753 }
754 }
755 }// triggered event
756
757 else { // not selected by physics selection task = not triggered
758 if(fUseMC && !fMcOnly){
759 if(CheckEvent(false)){
760
761 //read tracks
762 if(fESDEvent) ntracks = ReadEventESDMC(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 0);
763 else if(fAODEvent) ntracks = ReadEventAODMC(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 0);
28cc9dd5 764 else AliInfo("Fatal Error");
fa8e4941 765
3eb806f8 766 if (fCentralityMethod.Length() > 0)
11722f75 767 fNMcPrimAccept = centrality;
3eb806f8 768
fa8e4941 769 //analyse
770 if(pt.size()){
771 Analyse(pt, eta, phi, charge, strangenessWeight, fNMcPrimAccept, nTracksTracklets[1],nTracksTracklets[2], 0); //second part of step 0 // step 0 = AllAllMcNmc
772 }
773 }
774 }
775 }
776
777 if(fMcOnly){
778 // read event
779 if(fMode==0) ntracks = ReadEventESDMC(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 0);
780 else if (fMode==1) ntracks = ReadEventAODMC(pt, eta, phi, charge, strangenessWeight, nTracksTracklets, 0);
781
3eb806f8 782 if (fCentralityMethod.Length() > 0)
11722f75 783 fNMcPrimAccept = centrality;
3eb806f8 784
fa8e4941 785 // analyse
786 if(pt.size()){
787 Analyse(pt, eta, phi, charge, strangenessWeight, fNMcPrimAccept, nTracksTracklets[1],nTracksTracklets[2], 0);
788 }
8855d81c 789 }
117f99e3 790}
791
792
793//________________________________________________________________________
11722f75 794Double_t AliAnalysisTaskMinijet::ReadEventESD( vector<Float_t> &ptArray, vector<Float_t> &etaArray,
fa8e4941 795 vector<Float_t> &phiArray, vector<Short_t> &chargeArray,
796 vector<Float_t> &strangeArray,
11722f75 797 vector<Double_t> &nTracksTracklets, const Int_t step)
117f99e3 798{
fa8e4941 799 // gives back the number of esd tracks and pointer to arrays with track
800 // properties (pt, eta, phi)
801 // Uses TPC tracks with SPD vertex from now on
802
803 ptArray.clear();
804 etaArray.clear();
805 phiArray.clear();
806 chargeArray.clear();
807 strangeArray.clear();
808 nTracksTracklets.clear();
809
810 const AliESDVertex* vtxSPD = fESDEvent->GetPrimaryVertexSPD(); // uses track or SPD vertexer
811 fVertexZ[step]->Fill(vtxSPD->GetZ());
812
813 // Retreive the number of all tracks for this event.
814 Int_t ntracks = fESDEvent->GetNumberOfTracks();
815 if(fDebug>1) Printf("all ESD tracks: %d", ntracks);
7099c89a 816
fa8e4941 817 //first loop to check how many tracks are accepted
818 //------------------
11722f75 819 Double_t nAcceptedTracks=0;
fa8e4941 820 //Float_t nAcceptedTracksStrange=0;
821 for (Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
822
823 AliESDtrack *esdTrack = (AliESDtrack *)fESDEvent->GetTrack(iTracks);
824 if (!esdTrack) {
825 Error("ReadEventESD", "Could not receive track %d", iTracks);
826 continue;
827 }
828
829 // use TPC only tracks with non default SPD vertex
830 AliESDtrack *track = AliESDtrackCuts::
831 GetTPCOnlyTrack(dynamic_cast<AliESDEvent*>(fESDEvent),esdTrack->GetID());
832 if(!track) continue;
833 if(!fCuts->AcceptTrack(track)) {
834 delete track;
835 continue;// apply TPC track cuts before constrain to SPD vertex
836 }
837 if(track->Pt()>0.){
838 // only constrain tracks above threshold
839 AliExternalTrackParam exParam;
840 // take the B-field from the ESD, no 3D fieldMap available at this point
841 Bool_t relate = false;
842 relate = track->RelateToVertexTPC(vtxSPD,fESDEvent->GetMagneticField(),
843 kVeryBig,&exParam);
844 if(!relate){
845 delete track;
846 continue;
847 }
848 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),
849 exParam.GetCovariance());
850 }
851 else{
852 delete track;
853 continue;// only if tracks have pt<=0
854 }
855
856 if (TMath::Abs(track->Eta())<fEtaCut && track->Pt()>fPtMin && track->Pt()<fPtMax) {
857 ptArray.push_back(track->Pt());
858 etaArray.push_back(track->Eta());
859 phiArray.push_back(track->Phi());
860 chargeArray.push_back(track->Charge());
861 strangeArray.push_back(1);
862 ++nAcceptedTracks;
863 fHistPt->Fill(track->Pt());
864 }
865
866 // TPC only track memory needs to be cleaned up
867 if(track)
868 delete track;
869
7099c89a 870 }
fa8e4941 871
872 //need to be checked
873 if(nAcceptedTracks==0) return 0;
874
875 //tracklet loop
876 Int_t ntrackletsAccept=0;
877 AliMultiplicity * mult = (AliMultiplicity*)(fESDEvent->GetMultiplicity());
878 Int_t ntracklets = mult->GetNumberOfTracklets();
879 for(Int_t i=0;i< ntracklets;i++){
880 if(mult->GetDeltaPhi(i)<0.05 && TMath::Abs(mult->GetEta(i))<fEtaCut){
881 ntrackletsAccept++;
882 }
9ed461df 883 }
fa8e4941 884 nTracksTracklets.push_back(nAcceptedTracks);
885 nTracksTracklets.push_back(ntrackletsAccept);
886 nTracksTracklets.push_back(nAcceptedTracks);//in order to be compatible with mc analysis
887 //where here also neutral particles are counted.
888
889
890 fVzEvent=vtxSPD->GetZ(); // needed for correction map
891 if(step==5 || step ==2) {
892 fNRecAccept=nAcceptedTracks;
893 fNRecAcceptTracklet=ntrackletsAccept;
7099c89a 894 }
fa8e4941 895 return fNRecAccept;
896
897
898}
9ed461df 899
fa8e4941 900//________________________________________________________________________
11722f75 901Double_t AliAnalysisTaskMinijet::ReadEventESDRecMcProp( vector<Float_t> &ptArray, vector<Float_t> &etaArray,
fa8e4941 902 vector<Float_t> &phiArray, vector<Short_t> &chargeArray,
903 vector<Float_t> &strangeArray,
11722f75 904 vector<Double_t> &nTracksTracklets, const Int_t step)
fa8e4941 905{
906 // gives back the number of esd tracks and pointer to arrays with track
907 // properties (pt, eta, phi) of mc particles if available
908 // Uses TPC tracks with SPD vertex from now on
909
910 ptArray.clear();
911 etaArray.clear();
912 phiArray.clear();
913 chargeArray.clear();
914 strangeArray.clear();
915 nTracksTracklets.clear();
916
917
918 AliMCEvent *mcEvent = (AliMCEvent*) MCEvent();
919 if (!mcEvent) {
920 Error("ReadEventESDRecMcProp", "Could not retrieve MC event");
921 return 0;
7099c89a 922 }
fa8e4941 923 AliStack* stack = MCEvent()->Stack();
924 if(!stack) return 0;
9ed461df 925
fa8e4941 926
927 // Retreive the number of all tracks for this event.
928 Int_t ntracks = fESDEvent->GetNumberOfTracks();
929 if(fDebug>1)Printf("all ESD tracks: %d", ntracks);
930
931 const AliESDVertex *vtxSPD = fESDEvent->GetPrimaryVertexSPD();
932 fVertexZ[step]->Fill(vtxSPD->GetZ());
933
934 //track loop
11722f75 935 Double_t nAcceptedTracks=0;
fa8e4941 936 for (Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
937
938 AliVParticle *vtrack = fESDEvent->GetTrack(iTracks);
939 AliESDtrack *esdTrack = (AliESDtrack *)fESDEvent->GetTrack(iTracks);
940 if (!esdTrack) {
941 Error("ReadEventESDRecMcProp", "Could not receive track %d", iTracks);
942 continue;
943 }
944
945 // use TPC only tracks with non default SPD vertex
946 AliESDtrack *track = AliESDtrackCuts::
947 GetTPCOnlyTrack(dynamic_cast<AliESDEvent*>(fESDEvent),esdTrack->GetID());
948 if(!track) continue;
949 if(!fCuts->AcceptTrack(track)) {
950 delete track;
951 continue;// apply TPC track cuts before constrain to SPD vertex
952 }
953 if(track->Pt()>0.){
954 // only constrain tracks above threshold
955 AliExternalTrackParam exParam;
956 // take the B-field from the ESD, no 3D fieldMap available at this point
957 Bool_t relate = false;
958 relate = track->RelateToVertexTPC(vtxSPD,fESDEvent->GetMagneticField(),
959 kVeryBig,&exParam);
960 if(!relate){
961 delete track;
962 continue;
963 }
964 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),
965 exParam.GetCovariance());
966 }
967 else{
968 delete track;
969 continue;// only if tracks have pt<=0
970 }
971
972 //count tracks, if available, use mc particle properties
973 if(vtrack->GetLabel()<=0){
974 if (TMath::Abs(track->Eta())<fEtaCut && track->Pt()>fPtMin && track->Pt()<fPtMax){
975 ptArray.push_back(track->Pt());
976 etaArray.push_back(track->Eta());
977 phiArray.push_back(track->Phi());
978 chargeArray.push_back(track->Charge());
979 strangeArray.push_back(1);
980 ++nAcceptedTracks;
981 }
982 }
983 else{
984 TParticle *partOfTrack = stack->Particle(vtrack->GetLabel());
985 if (TMath::Abs(partOfTrack->Eta())<fEtaCut && partOfTrack->Pt()>fPtMin && partOfTrack->Pt()<fPtMax) {
986 ptArray.push_back(partOfTrack->Pt());
987 etaArray.push_back(partOfTrack->Eta());
988 phiArray.push_back(partOfTrack->Phi());
989 chargeArray.push_back(vtrack->Charge());
990 strangeArray.push_back(1);
991 ++nAcceptedTracks;
992 }
993 }
994
995 // TPC only track memory needs to be cleaned up
996 if(track)
997 delete track;
998
7099c89a 999 }
fa8e4941 1000
1001 if(nAcceptedTracks==0) return 0;
1002
1003 //tracklet loop
1004 Int_t ntrackletsAccept=0;
1005 AliMultiplicity * mult = (AliMultiplicity*)(fESDEvent->GetMultiplicity());
1006 Int_t ntracklets = mult->GetNumberOfTracklets();
1007 for(Int_t i=0;i< ntracklets;i++){
1008 if(mult->GetDeltaPhi(i)<0.05 && TMath::Abs(mult->GetEta(i))<fEtaCut){
1009 ntrackletsAccept++;
1010 }
1011 }
1012
1013 nTracksTracklets.push_back(nAcceptedTracks);
1014 nTracksTracklets.push_back(ntrackletsAccept);
1015 nTracksTracklets.push_back(nAcceptedTracks);//in order to be compatible with mc analysis
1016 //where here also neutral particles are counted.
1017
1018
1019 //get mc vertex for correction maps
1020 AliGenEventHeader* header = MCEvent()->GenEventHeader();
1021 TArrayF mcV;
1022 header->PrimaryVertex(mcV);
1023 fVzEvent= mcV[2];
1024
1025 return fNRecAccept; // give back reconstructed value
1026
1027
1028}
9ed461df 1029
117f99e3 1030
117f99e3 1031
9ed461df 1032
7099c89a 1033//________________________________________________________________________
11722f75 1034Double_t AliAnalysisTaskMinijet::ReadEventESDMC(vector<Float_t> &ptArray, vector<Float_t> &etaArray,
fa8e4941 1035 vector<Float_t> &phiArray, vector<Short_t> &chargeArray,
1036 vector<Float_t> &strangeArray,
11722f75 1037 vector<Double_t> &nTracksTracklets, const Int_t step)
fa8e4941 1038{
1039 // gives back the number of charged prim MC particle and pointer to arrays
1040 // with particle properties (pt, eta, phi)
1041
1042 ptArray.clear();
1043 etaArray.clear();
1044 phiArray.clear();
1045 chargeArray.clear();
1046 strangeArray.clear();
1047 nTracksTracklets.clear();
1048
1049 fNMcPrimAccept=0;
1050
1051 AliMCEvent *mcEvent = (AliMCEvent*) MCEvent();
1052 if (!mcEvent) {
1053 Error("ReadEventESDMC", "Could not retrieve MC event");
1054 return 0;
117f99e3 1055 }
7099c89a 1056
fa8e4941 1057 AliStack* stack = MCEvent()->Stack();
1058 if(!stack) return 0;
1059
1060 Int_t ntracks = mcEvent->GetNumberOfTracks();
1061 if(fDebug>1)Printf("MC particles: %d", ntracks);
1062
1063 //vertex
1064 AliGenEventHeader* header = MCEvent()->GenEventHeader();
1065 TArrayF mcV;
1066 Float_t vzMC=0.;
1067 if(header){
1068 header->PrimaryVertex(mcV);
1069 vzMC = mcV[2];
1070 if(step==1){
1071 fVertexZ[0]->Fill(vzMC);
1072 }
1073 fVertexZ[step]->Fill(vzMC);
7099c89a 1074 }
fa8e4941 1075
1076 //----------------------------------
1077 //first track loop to check how many chared primary tracks are available
11722f75 1078 Double_t nChargedPrimaries=0;
1079 Double_t nAllPrimaries=0;
fa8e4941 1080
11722f75 1081 Double_t nPseudoTracklets=0;
fa8e4941 1082 for (Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
1083 AliMCParticle *track = dynamic_cast<AliMCParticle*>(mcEvent->GetTrack(iTracks));
1084 if (!track) {
1085 Error("ReadEventESDMC", "Could not receive track %d", iTracks);
1086 continue;
1087 }
1088
1089
1090 if(//count also charged particles in case of fSelectParticles==2 (only neutral)
1091 !SelectParticlePlusCharged(
1092 track->Charge(),
1093 track->PdgCode(),
1094 stack->IsPhysicalPrimary(track->Label())
1095 )
1096 )
1097 continue;
1098
1099 //count number of pseudo tracklets
1100 if(TMath::Abs(track->Eta())<=fEtaCut && track->Pt()>0.0) ++nPseudoTracklets; //0.035
1101 //same cuts as on ESDtracks
1102 if(TMath::Abs(track->Eta())>fEtaCut || track->Pt()<fPtMin || track->Pt()>fPtMax) continue;
1103
1104 //count all primaries
1105 ++nAllPrimaries;
1106 //count charged primaries
1107 if (track->Charge() != 0) ++nChargedPrimaries;
1108
1109 if(fDebug>2) Printf("PDG=%d, IsPrim=%d", track->PdgCode(),stack->IsPhysicalPrimary(track->Label()) );
1110
1111
9ed461df 1112 }
fa8e4941 1113 //----------------------------------
1114
1115 if(fDebug>2){
11722f75 1116 Printf("All in acceptance=%f",nAllPrimaries);
1117 Printf("Charged in acceptance =%f",nChargedPrimaries);
7099c89a 1118 }
fa8e4941 1119
1120 fChargedPi0->Fill(nAllPrimaries,nChargedPrimaries);
1121
1122 if(nAllPrimaries==0) return 0;
1123 if(nChargedPrimaries==0) return 0;
1124
1125 //track loop
1126 //Int_t iChargedPiK=0;
1127 for (Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
1128 AliMCParticle *track = dynamic_cast<AliMCParticle*>(mcEvent->GetTrack(iTracks));
1129 if (!track) {
1130 Error("ReadEventESDMC", "Could not receive track %d", iTracks);
1131 continue;
1132 }
1133
1134 if(!SelectParticle(
1135 track->Charge(),
1136 track->PdgCode(),
1137 stack->IsPhysicalPrimary(track->Label())
1138 )
1139 )
1140 continue;
1141
1142
1143 //same cuts as on ESDtracks
1144 if(TMath::Abs(track->Eta())>fEtaCut || track->Pt()<fPtMin || track->Pt()>fPtMax) continue;
1145
1146 if(fDebug>2) Printf("After: PDG=%d, IsPrim=%d", track->PdgCode(),stack->IsPhysicalPrimary(track->Label()) );
1147
1148
1149 fHistPtMC->Fill(track->Pt());
1150 //fills arrays with track properties
1151 ptArray.push_back(track->Pt());
1152 etaArray.push_back(track->Eta());
1153 phiArray.push_back(track->Phi());
1154 chargeArray.push_back(track->Charge());
1155 strangeArray.push_back(1);
1156
1157 } //track loop
1158
1159 nTracksTracklets.push_back(nChargedPrimaries);
1160 nTracksTracklets.push_back(nPseudoTracklets);
1161 if(fSelectParticles!=2){
1162 nTracksTracklets.push_back(nAllPrimaries);
117f99e3 1163 }
7099c89a 1164 else{
fa8e4941 1165 nTracksTracklets.push_back(nAllPrimaries-nChargedPrimaries); // only neutral
7099c89a 1166 }
fa8e4941 1167
1168 fNMcPrimAccept = nChargedPrimaries;
1169 fNMcPrimAcceptTracklet = nPseudoTracklets;
1170
9cccc98a 1171 if(step==1){
fa8e4941 1172 fNmcNch->Fill(fNMcPrimAccept,fNRecAccept);
1173 fPNmcNch->Fill(fNMcPrimAccept,fNRecAccept);
1174 fNmcNchTracklet->Fill(fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
1175 fPNmcNchTracklet->Fill(fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
9cccc98a 1176 }
fa8e4941 1177 if(step==3){
1178 fNmcNchVtx->Fill(fNMcPrimAccept,fNRecAccept);
1179 fPNmcNchVtx->Fill(fNMcPrimAccept,fNRecAccept);
1180 fNmcNchVtxTracklet->Fill(fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
1181 fPNmcNchVtxTracklet->Fill(fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
117f99e3 1182 }
117f99e3 1183
fa8e4941 1184 fVzEvent= vzMC;
1185 return fNRecAccept;
7099c89a 1186
117f99e3 1187}
1188
1189//________________________________________________________________________
11722f75 1190Double_t AliAnalysisTaskMinijet::ReadEventAOD( vector<Float_t> &ptArray, vector<Float_t> &etaArray,
fa8e4941 1191 vector<Float_t> &phiArray, vector<Short_t> &chargeArray,
1192 vector<Float_t> &strangeArray,
11722f75 1193 vector<Double_t> &nTracksTracklets, const Int_t step)
117f99e3 1194{
fa8e4941 1195 // gives back the number of AOD tracks and pointer to arrays with track
1196 // properties (pt, eta, phi)
1197
1198 ptArray.clear();
1199 etaArray.clear();
1200 phiArray.clear();
1201 chargeArray.clear();
1202 strangeArray.clear();
1203 nTracksTracklets.clear();
1204
1205 TClonesArray *mcArray=0x0;
1206 if(fAnalysePrimOnly || (fCorrStrangeness && fUseMC)){
1207 mcArray = (TClonesArray*)fAODEvent->FindListObject(AliAODMCParticle::StdBranchName());
117f99e3 1208 }
fa8e4941 1209
1210
1211 AliAODVertex* vertex= (AliAODVertex*)fAODEvent->GetPrimaryVertexSPD();//GetPrimaryVertex()
1212 Double_t vzAOD=vertex->GetZ();
1213 fVertexZ[step]->Fill(vzAOD);
1214
1215 // Retreive the number of tracks for this event.
1216 Int_t ntracks = fAODEvent->GetNumberOfTracks();
1217 if(fDebug>1) Printf("AOD tracks: %d", ntracks);
1218
1219
11722f75 1220 Double_t nAcceptedTracks=0;
fa8e4941 1221 Float_t nAcceptedTracksStrange=0;
1222 for (Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
1223 AliAODTrack *track = (AliAODTrack *)fAODEvent->GetTrack(iTracks);
1224 if (!track) {
1225 Error("ReadEventAOD", "Could not receive track %d", iTracks);
1226 continue;
1227 }
1228
1229 AliVParticle *vtrack = fAODEvent->GetTrack(iTracks);
1230
1231 //use only tracks from primaries
1232 if(fAnalysePrimOnly){
1233 if(vtrack->GetLabel()<=0)continue;
1234 if(!(static_cast<AliAODMCParticle*>(mcArray->At(vtrack->GetLabel()))->IsPhysicalPrimary()))continue;
1235 }
1236
1237 Double_t save= track->Pt();
1238 Double_t d0rphiz[2],covd0[3];
14fb8e86 1239
1240 AliAODTrack* clone = (AliAODTrack*) track->Clone("trk_clone"); //need clone, in order not to change track parameters
1241 Bool_t isDca= clone->PropagateToDCA(fAODEvent->GetPrimaryVertex(),fAODEvent->GetMagneticField(),9999.,d0rphiz,covd0);
1242 delete clone;
fa8e4941 1243 fPropagateDca->Fill(Int_t(isDca));
1244 if(TMath::Abs(save - track->Pt())>1e-6) Printf("Before pt=%f, After pt=%f",save, track->Pt());
1245
1246 if(track->TestFilterBit(fFilterBit) && TMath::Abs(track->Eta())<fEtaCut
1247 && track->Pt()>fPtMin && track->Pt()<fPtMax){
1248
1249 nAcceptedTracks++;
1250
1251 // Printf("dca= %f", track->DCA());
1252 //save track properties in vector
1253 ptArray.push_back(track->Pt());
1254 etaArray.push_back(track->Eta());
1255 phiArray.push_back(track->Phi());
1256 chargeArray.push_back(track->Charge());
1257 fHistPt->Fill(track->Pt());
1258
1259
1260 //correction for underestimation of strangeness in Monte Carlos -> underestimation of contamination
1261 Float_t factor=1.;
1262 if(fUseMC && fCorrStrangeness && step==7){
1263 if(vtrack->GetLabel()>0){
1264 AliAODMCParticle* mcprong =(AliAODMCParticle*)mcArray->At(vtrack->GetLabel());
1265 if(mcprong){
1266 Int_t labmom = mcprong->GetMother();
1267 if(labmom>=0){
1268 AliAODMCParticle* mcmother=(AliAODMCParticle*)mcArray->At(labmom);
1269 Int_t pdgMother=0;
1270 if(mcmother) {
1271 pdgMother = mcmother->GetPdgCode();
1272 if(TMath::Abs(pdgMother)==310 || TMath::Abs(pdgMother)==130 || TMath::Abs(pdgMother)==321){ //K^0_S, K^0_L, K^+-
1273 if(track->Pt()<=1)factor=1./0.7; // values from strangeness publication
1274 else factor=1./0.6;// values from strangeness publication
1275 }
1276 if(TMath::Abs(pdgMother)==3122) { //Lambda
1277 factor=1./0.25; // values from strangeness publication
1278 }
1279 }
1280 }
1281 }
1282 }
1283 }
1284 nAcceptedTracksStrange+=factor;
1285 strangeArray.push_back(factor);
1286 fDcaXY[step]->Fill(d0rphiz[0], factor);
1287 fDcaZ[step]->Fill(d0rphiz[0], factor);
1288
1289 }
7099c89a 1290 }
fa8e4941 1291 //need to check this option for MC
1292 if(nAcceptedTracks==0) return 0;
7099c89a 1293
2b52b282 1294
fa8e4941 1295 //tracklet loop
1296 Int_t ntrackletsAccept=0;
1297 AliAODTracklets * mult= (AliAODTracklets*)fAODEvent->GetTracklets();
1298 for(Int_t i=0;i<mult->GetNumberOfTracklets();++i){
1299 if(TMath::Abs(mult->GetDeltaPhi(i))<0.05 && TMath::Abs(TMath::Log(TMath::Tan(0.5 * mult->GetTheta(i))))<fEtaCut){
1300 ++ntrackletsAccept;
1301 }
7099c89a 1302 }
fa8e4941 1303
1304
1305 nTracksTracklets.push_back(nAcceptedTracks);
1306 nTracksTracklets.push_back(ntrackletsAccept);
1307 nTracksTracklets.push_back(nAcceptedTracks);//in order to be compatible with mc analysis
1308 //where here also neutral particles are counted.
1309
1310
1311 fVzEvent= vzAOD;
1312 if(step==5 || step==2){
1313 fNRecAccept = nAcceptedTracks; // needed for MC case //step5 = TrigVtxRecNrec
1314 fNRecAcceptTracklet = ntrackletsAccept; // needed for MC case //step5 = TrigVtxRecNrec
7099c89a 1315 }
fa8e4941 1316 if(step==7)fNRecAcceptStrangeCorr = nAcceptedTracksStrange;
1317
1318 return fNRecAccept; // at the moment, always return reconstructed multiplicity
1319
7099c89a 1320}
1321
1322//________________________________________________________________________
11722f75 1323Double_t AliAnalysisTaskMinijet::ReadEventAODRecMcProp( vector<Float_t> &ptArray, vector<Float_t> &etaArray,
fa8e4941 1324 vector<Float_t> &phiArray, vector<Short_t> &chargeArray,
1325 vector<Float_t> &strangeArray,
11722f75 1326 vector<Double_t> &nTracksTracklets, const Int_t step)
7099c89a 1327{
fa8e4941 1328 // gives back the number of AOD tracks and pointer to arrays with track
1329 // properties (pt, eta, phi)
1330
1331
1332 ptArray.clear();
1333 etaArray.clear();
1334 phiArray.clear();
1335 chargeArray.clear();
1336 strangeArray.clear();
1337 nTracksTracklets.clear();
1338
1339
1340 // Retreive the number of tracks for this event.
1341 Int_t ntracks = fAODEvent->GetNumberOfTracks();
1342 if(fDebug>1) Printf("AOD tracks: %d", ntracks);
1343
1344
1345 //get array of mc particles
1346 TClonesArray *mcArray = (TClonesArray*)fAODEvent->
7099c89a 1347 FindListObject(AliAODMCParticle::StdBranchName());
fa8e4941 1348 if(!mcArray){
28cc9dd5 1349 AliInfo("No MC particle branch found");
fa8e4941 1350 return kFALSE;
117f99e3 1351 }
fa8e4941 1352
1353 AliAODVertex* vtx= (AliAODVertex*)fAODEvent->GetPrimaryVertexSPD();//GetPrimaryVertex()
1354 Double_t vzAOD=vtx->GetZ();
1355 fVertexZ[step]->Fill(vzAOD);
1356
11722f75 1357 Double_t nAcceptedTracks=0;
fa8e4941 1358 for (Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
1359 AliAODTrack *track = (AliAODTrack *)fAODEvent->GetTrack(iTracks);
1360
1361 AliVParticle *vtrack = fAODEvent->GetTrack(iTracks);
1362
1363 if (!track) {
1364 Error("ReadEventAODRecMcProp", "Could not receive track %d", iTracks);
1365 continue;
1366 }
1367
1368 //use only tracks from primaries
1369 if(fAnalysePrimOnly){
1370 if(vtrack->GetLabel()<=0)continue;
1371 if(!(static_cast<AliAODMCParticle*>(mcArray->At(vtrack->GetLabel()))->IsPhysicalPrimary()))continue;
1372 }
1373
1374 if(track->TestFilterBit(fFilterBit) && TMath::Abs(track->Eta())<fEtaCut &&
1375 track->Pt()>fPtMin && track->Pt()<fPtMax){
1376
1377 nAcceptedTracks++;
1378 Float_t factor=1.;
1379
1380 //save track properties in vector
1381 if(vtrack->GetLabel()<=0){ //fake tracks before "label<0", but crash in AOD079 // what is the meaning of label 0
1382 // Printf("Fake track");
1383 // continue;
1384 ptArray.push_back(track->Pt());
1385 etaArray.push_back(track->Eta());
1386 phiArray.push_back(track->Phi());
1387 chargeArray.push_back(track->Charge());
2b52b282 1388
fa8e4941 1389 }
1390 else{//mc properties
1391 AliAODMCParticle *partOfTrack = (AliAODMCParticle*)mcArray->At(vtrack->GetLabel());
1392 if(!partOfTrack) Printf("label=%d", vtrack->GetLabel());
1393 if(partOfTrack){
1394 ptArray.push_back(partOfTrack->Pt());
1395 etaArray.push_back(partOfTrack->Eta());
1396 phiArray.push_back(partOfTrack->Phi());
1397 chargeArray.push_back(vtrack->Charge());//partOfTrack?
1398
1399 //correction for underestimation of strangeness in Monte Carlos -> underestimation of contamination
1400 if(fUseMC && fCorrStrangeness && step==6){
1401 Int_t labmom = partOfTrack->GetMother();
1402 if(labmom>=0){
1403 AliAODMCParticle* mcmother=(AliAODMCParticle*)mcArray->At(labmom);
1404 Int_t pdgMother=0;
1405 if(mcmother) {
1406 pdgMother = mcmother->GetPdgCode();
1407 if(TMath::Abs(pdgMother)==310 || TMath::Abs(pdgMother)==130 || TMath::Abs(pdgMother)==321){ //K^0_S, K^0_L, K^+-
1408 if(track->Pt()<=1)factor=1./0.7; // values from strangeness publication
1409 else factor=1./0.6;// values from strangeness publication
1410 }
1411 if(TMath::Abs(pdgMother)==3122) { //Lambda
1412 factor=1./0.25; // values from strangeness publication
1413 }
1414 }
1415 }
1416 }
1417 }
1418 }
1419 strangeArray.push_back(factor);
1420
1421 }
7099c89a 1422 }
fa8e4941 1423 //need to check this option for MC
1424 if(nAcceptedTracks==0) return 0;
1425
1426 //tracklet loop
1427 Int_t ntrackletsAccept=0;
1428 AliAODTracklets * mult= (AliAODTracklets*)fAODEvent->GetTracklets();
1429 for(Int_t i=0;i<mult->GetNumberOfTracklets();++i){
1430 if(TMath::Abs(mult->GetDeltaPhi(i))<0.05 && TMath::Abs(TMath::Log(TMath::Tan(0.5 * mult->GetTheta(i))))<fEtaCut ){
1431 ++ntrackletsAccept;
1432 }
117f99e3 1433 }
fa8e4941 1434
1435
1436 nTracksTracklets.push_back(nAcceptedTracks);
1437 nTracksTracklets.push_back(ntrackletsAccept);
1438 nTracksTracklets.push_back(nAcceptedTracks);//in order to be compatible with mc analysis
1439 //where here also neutral particles are counted.
1440
1441
1442 //check vertex (mc)
1443 AliAODMCHeader *aodMCheader = (AliAODMCHeader *) fAODEvent->
7099c89a 1444 FindListObject(AliAODMCHeader::StdBranchName());
fa8e4941 1445 Float_t vzMC = aodMCheader->GetVtxZ();
1446
1447 fVzEvent= vzMC;
1448 return fNRecAccept;//this is the rec value from step 5
1449
1450}
7099c89a 1451
1452
117f99e3 1453
1454//________________________________________________________________________
11722f75 1455Double_t AliAnalysisTaskMinijet::ReadEventAODMC( vector<Float_t> &ptArray, vector<Float_t> &etaArray,
fa8e4941 1456 vector<Float_t> &phiArray, vector<Short_t> &chargeArray,
1457 vector<Float_t> &strangeArray,
11722f75 1458 vector<Double_t> &nTracksTracklets, const Int_t step)
117f99e3 1459{
fa8e4941 1460 // gives back the number of AOD MC particles and pointer to arrays with particle
1461 // properties (pt, eta, phi)
1462
1463 ptArray.clear();
1464 etaArray.clear();
1465 phiArray.clear();
1466 chargeArray.clear();
1467 strangeArray.clear();
1468 nTracksTracklets.clear();
1469
1470 //check vertex
1471 AliAODMCHeader *aodMCheader = (AliAODMCHeader *) fAODEvent->
7099c89a 1472 FindListObject(AliAODMCHeader::StdBranchName());
fa8e4941 1473 Float_t vzMC = aodMCheader->GetVtxZ();
1474 if(step==1){
1475 fVertexZ[0]->Fill(vzMC);
1476 }
1477 fVertexZ[step]->Fill(vzMC);
1478
1479 //retreive MC particles from event
1480 TClonesArray *mcArray = (TClonesArray*)fAODEvent->
b7c0438e 1481 FindListObject(AliAODMCParticle::StdBranchName());
fa8e4941 1482 if(!mcArray){
28cc9dd5 1483 AliInfo("No MC particle branch found");
fa8e4941 1484 return kFALSE;
117f99e3 1485 }
fa8e4941 1486
1487 Int_t ntracks = mcArray->GetEntriesFast();
1488 if(fDebug>1)Printf("MC particles: %d", ntracks);
1489
1490
1491 // Track loop: chek how many particles will be accepted
1492 //Float_t vzMC=0.;
11722f75 1493 Double_t nChargedPrim=0;
1494 Double_t nAllPrim=0;
1495 Double_t nPseudoTracklets=0;
fa8e4941 1496 for (Int_t it = 0; it < ntracks; it++) {
1497 AliAODMCParticle *track = (AliAODMCParticle*)mcArray->At(it);
1498 if (!track) {
1499 Error("ReadEventAODMC", "Could not receive particle %d", it);
1500 continue;
1501 }
1502
1503 if(!SelectParticlePlusCharged(
1504 track->Charge(),
1505 track->PdgCode(),
1506 track->IsPhysicalPrimary()
1507 )
1508 )
1509 continue;
1510
1511 if(TMath::Abs(track->Eta())<fEtaCut && track->Pt()>0.0)++nPseudoTracklets; //0.035
1512 if(TMath::Abs(track->Eta())>fEtaCut || track->Pt()<fPtMin || track->Pt()>fPtMax) continue;
1513
1514 nAllPrim++;
1515 if(track->Charge()!=0) nChargedPrim++;
1516
b7c0438e 1517 }
fa8e4941 1518
1519
1520 if(nAllPrim==0) return 0;
1521 if(nChargedPrim==0) return 0;
1522
1523 //generate array with size of number of accepted tracks
1524 fChargedPi0->Fill(nAllPrim,nChargedPrim);
1525
1526
1527 // Track loop: fill arrays for accepted tracks
1528 // Int_t iChargedPrim=0;
1529 for (Int_t it = 0; it < ntracks; it++) {
1530 AliAODMCParticle *track = (AliAODMCParticle*)mcArray->At(it);
1531 if (!track) {
1532 Error("ReadEventAODMC", "Could not receive particle %d", it);
1533 continue;
1534 }
1535
1536 if(!SelectParticle(
1537 track->Charge(),
1538 track->PdgCode(),
1539 track->IsPhysicalPrimary()
1540 )
1541 )
1542 continue;
1543
1544 if(TMath::Abs(track->Eta())>fEtaCut || track->Pt()<fPtMin || track->Pt()>fPtMax) continue;
1545
1546 fHistPtMC->Fill(track->Pt());
1547 ptArray.push_back(track->Pt());
1548 etaArray.push_back(track->Eta());
1549 phiArray.push_back(track->Phi());
1550 chargeArray.push_back(track->Charge());
1551 strangeArray.push_back(1);
1552 }
1553
1554 nTracksTracklets.push_back(nChargedPrim);
1555 nTracksTracklets.push_back(nPseudoTracklets);
1556 if(fSelectParticles!=2){
1557 nTracksTracklets.push_back(nAllPrim);
1558 }
1559 else{
1560 nTracksTracklets.push_back(nAllPrim-nChargedPrim); // only neutral
1561 }
1562
1563
1564
1565 fVzEvent= vzMC;
1566 fNMcPrimAccept = nChargedPrim;
1567 fNMcPrimAcceptTracklet = nPseudoTracklets;
1568
1569 if(step==1){ // step 1 = Trig All Mc Nmc
1570 fNmcNch->Fill( fNMcPrimAccept,fNRecAccept);
1571 fPNmcNch->Fill(fNMcPrimAccept,fNRecAccept);
1572 fNmcNchTracklet->Fill( fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
1573 fPNmcNchTracklet->Fill(fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
1574 }
1575 if(step==3){ // step 3 = Trig vtx Mc
1576 fNmcNchVtx->Fill( fNMcPrimAccept,fNRecAccept);
1577 fNmcNchVtxStrangeCorr->Fill( fNMcPrimAccept,fNRecAcceptStrangeCorr);
1578 fPNmcNchVtx->Fill(fNMcPrimAccept,fNRecAccept);
1579 fNmcNchVtxTracklet->Fill( fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
1580 fPNmcNchVtxTracklet->Fill(fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
1581 }
1582 return fNRecAccept; // rec value from step 5 or step 2
1583
1584
117f99e3 1585}
1586
1587//________________________________________________________________________
fa8e4941 1588void AliAnalysisTaskMinijet::Analyse(const vector<Float_t> &pt,
1589 const vector<Float_t> &eta,
1590 const vector<Float_t> &phi,
1591 const vector<Short_t> &charge,
1592 const vector<Float_t> &strangeWeight,
11722f75 1593 const Double_t ntracksCharged,
fa8e4941 1594 const Int_t ntracklets,
1595 const Int_t nAll,
1596 const Int_t step)
117f99e3 1597{
fa8e4941 1598
1599 // analyse track properties (comming from either ESDs or AODs) in order to compute
1600 // mini jet activity (chared tracks) as function of charged multiplicity
1601
1602 fStep->Fill(step);
1603
1604 if(fDebug){
1605 Printf("Analysis Step=%d", step);
1606 if(fDebug>2){
1607 Printf("nAll=%d",nAll);
11722f75 1608 Printf("nCharged=%f",ntracksCharged);
fa8e4941 1609 for (Int_t i = 0; i < nAll; i++) {
1610 Printf("pt[%d]=%f",i,pt[i]);
1611 }
1612 }
7099c89a 1613 }
fa8e4941 1614
1615 //calculation of mean pt for all tracks in case of step==0
1616 if(step==5 || step ==2){ // rec step
1617 Double_t meanPt=0.;
1618 Double_t leadingPt=0.;
1619 for (Int_t i = 0; i < nAll; i++) {
1620 if(pt[i]<0.01)continue;
1621 meanPt+=pt[i];
1622 if(leadingPt<pt[i])leadingPt=pt[i];
1623 }
1624 meanPt=meanPt/nAll;
1625 fMeanPtRec=meanPt;
1626 fLeadingPtRec=leadingPt;
1627 }
1628
1629 Double_t propEvent[4] = {ntracksCharged,fVzEvent,fMeanPtRec, fLeadingPtRec}; //vz: {rec, mc, mc}, meanPt and Nrec is always rec value
1630 fMapEvent[step]->Fill(propEvent);
1631
1632 fNcharge[step]->Fill(ntracksCharged);
1633
1634 Float_t ptEventAxis=0; // pt event axis
1635 Float_t etaEventAxis=0; // eta event axis
1636 Float_t phiEventAxis=0; // phi event axis
1637 Short_t chargeEventAxis=0; // charge event axis
1638 Float_t strangeWeightEventAxis=0; // strange weight event axis
1639
1640 Float_t ptOthers = 0; // pt others // for all other tracks around event axis -> see loop
1641 Float_t etaOthers = 0; // eta others
1642 Float_t phiOthers = 0; // phi others
1643 Short_t chargeOthers = 0; // charge others
1644 Float_t strangeWeightOthers = 0; // strange weight others
1645
1646 Int_t *pindexInnerEta = new Int_t [nAll+1];
1647 Float_t *ptInnerEta = new Float_t[nAll+1];
1648
1649
1650
7099c89a 1651 for (Int_t i = 0; i < nAll; i++) {
fa8e4941 1652
1653 if(pt[i]<0.01)continue;
1654
1655 //fill single particle correction for first step of pair correction
1656 Double_t propAll[3] = {eta[i],pt[i],ntracksCharged };
1657 fMapAll[step]->Fill(propAll, strangeWeight[i]);
1658
1659
1660 //filling of simple check plots
1661 if(pt[i]<0.7) continue;
1662 fPt[step] -> Fill( pt[i]);
1663 fEta[step] -> Fill(eta[i]);
1664 fPhi[step] -> Fill(phi[i]);
1665 fPhiEta[step]-> Fill(phi[i], eta[i]);
1666
1667 pindexInnerEta[i]=0; //set all values to zero
1668 //fill new array for eta check
1669 ptInnerEta[i]= pt[i];
1670
7099c89a 1671 }
fa8e4941 1672
1673
1674
1675 // define event axis: leading or random track (with pt>fTriggerPtCut)
1676 // ---------------------------------------
1677 Int_t highPtTracks=0;
1678 Int_t highPtTracksInnerEta=0;
1679 // Double_t highPtTracksInnerEtaCorr=0;
1680 Int_t mult09=0;
1681
1682 //count high pt tracks and high pt tracks in acceptance for seeds
1683 for(Int_t i=0;i<nAll;i++){
1684
1685 if(pt[i]<0.01)continue;
1686
1687 if(TMath::Abs(eta[i])<0.9){
1688 mult09++;
1689 }
1690
1691 if(pt[i]>fTriggerPtCut) {
1692 highPtTracks++;
1693 }
1694
1695 // seed should be place in middle of acceptance, that complete cone is in acceptance
1696 if(pt[i]>fTriggerPtCut && TMath::Abs(eta[i])<fEtaCutSeed && charge[i]!=0){
1697
1698 // Printf("eta=%f", eta[i]);
1699 highPtTracksInnerEta++;
1700
1701 }
1702 else{
1703 ptInnerEta[i]=0;
1704 }
9ed461df 1705 }
fa8e4941 1706
1707
1708 //sort array in order to get highest pt tracks first
1709 //index can be computed with pindexInnerEta[number]
1710 if(nAll) TMath::Sort(nAll, ptInnerEta, pindexInnerEta, kTRUE);
1711
1712 // plot of multiplicity distributions
1713 fNch07Nch[step]->Fill(ntracksCharged, highPtTracksInnerEta);
1714 fPNch07Nch[step]->Fill(ntracksCharged, highPtTracksInnerEta);
1715
1716 if(ntracklets){
1717 fNch07Tracklet[step]->Fill(ntracklets, highPtTracksInnerEta);//only counts tracks which can be used as seeds
1718 fNchTracklet[step]->Fill(ntracklets, ntracksCharged);
1719 fPNch07Tracklet[step]->Fill(ntracklets, highPtTracksInnerEta);//only counts tracks which can be used as seeds
b9ad6f04 1720 }
fa8e4941 1721
1722 //analysis can only be performed with event axis, defined by high pt track
1723
1724
1725 if(highPtTracks>0 && highPtTracksInnerEta>0){
1726
1727 // build pairs in two track loops
1728 // loop over all possible trigger particles (defined by pt_trig and eta_acceptance)
1729 for(Int_t axis=0;(axis<nAll) && (pt[pindexInnerEta[axis]]>=fTriggerPtCut); axis++){
1730
1731 //EventAxisRandom track properties
1732 ptEventAxis = pt [pindexInnerEta[axis]];
1733 etaEventAxis = eta[pindexInnerEta[axis]];
1734 phiEventAxis = phi[pindexInnerEta[axis]];
1735 chargeEventAxis = charge[pindexInnerEta[axis]];
1736 strangeWeightEventAxis = strangeWeight[pindexInnerEta[axis]];
1737 fPtSeed[step] -> Fill( ptEventAxis);
1738 fEtaSeed[step] -> Fill(etaEventAxis);
1739 fPhiSeed[step] -> Fill(phiEventAxis);
1740
1741
1742 Double_t prop[3] = {etaEventAxis,ptEventAxis,ntracksCharged };
1743 fMapSingleTrig[step]->Fill(prop, strangeWeightEventAxis);
1744
1745 //associated tracks
1746 for (Int_t iTrack = axis+1; iTrack < nAll; iTrack++) {
1747
1748 if(pt[pindexInnerEta[iTrack]]<fAssociatePtCut) continue;
1749
1750 if(fSelectParticlesAssoc==1){
1751 if(charge[pindexInnerEta[iTrack]]==0)continue;
1752 }
1753 if(fSelectParticlesAssoc==2){
1754 if(charge[pindexInnerEta[iTrack]]!=0)continue;
1755 }
1756
1757
1758 ptOthers = pt [pindexInnerEta[iTrack]];
1759 etaOthers = eta[pindexInnerEta[iTrack]];
1760 phiOthers = phi[pindexInnerEta[iTrack]];
1761 chargeOthers = charge[pindexInnerEta[iTrack]];
1762 strangeWeightOthers = strangeWeight[pindexInnerEta[iTrack]];
1763
1764
1765 //plot only properties of tracks with pt>ptassoc
1766 fPtOthers[step] -> Fill( ptOthers);
1767 fEtaOthers[step] -> Fill(etaOthers);
1768 fPhiOthers[step] -> Fill(phiOthers);
1769 fPtEtaOthers[step] -> Fill(ptOthers, etaOthers);
1770
1771 // if(fDebug>2)Printf("%f, %f", pt[pindexInnerEta[axis]], pt[pindexInnerEta[iTrack]]);
1772
1773 Float_t dPhi = phiOthers-phiEventAxis;
1774 if(dPhi>1.5*TMath::Pi()) dPhi = dPhi-2*TMath::Pi();
1775 else if(dPhi<-0.5*TMath::Pi())dPhi=dPhi+2*TMath::Pi();
1776 Float_t dEta=etaOthers-etaEventAxis;
1777
1778
1779 fDPhiDEtaEventAxis[step]->Fill(dPhi, dEta, strangeWeightEventAxis*strangeWeightOthers);
1780 fDPhiEventAxis[step]->Fill(dPhi, strangeWeightEventAxis*strangeWeightOthers);
1781
1782 //check outliers
28cc9dd5 1783 if(ptEventAxis< 0.4 || ptEventAxis > 100) AliInfo("particles out of range pt");
1784 if(ntracksCharged<-1 || ntracksCharged>1500) AliInfo("particles out of range ncharge");
fa8e4941 1785 if(TMath::Abs(dEta)>2*fEtaCut) {
28cc9dd5 1786 AliInfo("particles out of range dEta");
c8ade419 1787 AliInfo(Form("eta1=%f, eta2=%f", etaOthers, etaEventAxis));
1788 AliInfo(Form("step=%d",step));
fa8e4941 1789 }
1790 if(dPhi<-0.5*TMath::Pi() || dPhi>1.5*TMath::Pi()){
c8ade419 1791 AliInfo(Form("particles out of range dPhi"));
1792 AliInfo(Form("phi1=%f, phi2=%f", phiOthers, phiEventAxis));
fa8e4941 1793 }
1794
1795 Bool_t isLikeSign = CheckLikeSign(chargeEventAxis, chargeOthers);
1796
2942f542 1797 Double_t prop6[6] = {ptEventAxis,ptOthers,dEta,dPhi,ntracksCharged, static_cast<Double_t>(isLikeSign) };
fa8e4941 1798 fMapPair[step]->Fill(prop6, strangeWeightEventAxis*strangeWeightOthers);
1799
1800 //thrid track loop (Andreas: three particle auto-correlations)
1801 if(fThreeParticleCorr){
1802 for (Int_t third = iTrack+1; third < nAll; third++) {
1803 if(pt[pindexInnerEta[iTrack]]<fTriggerPtCut) continue;
1804 if(pt[pindexInnerEta[third]]<fTriggerPtCut) continue;
1805
1806 //dphi1
1807 Float_t dPhi1 = phiEventAxis - phiOthers;
1808 if(dPhi1>1.5*TMath::Pi()) dPhi1 = dPhi1-2*TMath::Pi();
1809 else if(dPhi1<-0.5*TMath::Pi())dPhi1=dPhi1+2*TMath::Pi();
1810
1811 Float_t phiThird = phi[pindexInnerEta[third]];
1812 Float_t strangeWeightThird = strangeWeight[pindexInnerEta[third]];
1813
1814 //dphi2
1815 Float_t dPhi2 = phiEventAxis - phiThird;
1816 if(dPhi2>1.5*TMath::Pi()) dPhi2 = dPhi2-2*TMath::Pi();
1817 else if(dPhi2<-0.5*TMath::Pi())dPhi2=dPhi2+2*TMath::Pi();
1818
1819 fDPhi1DPhi2[step]-> Fill(dPhi1, dPhi2, strangeWeightEventAxis*strangeWeightOthers*strangeWeightThird);
1820 Double_t propThree[3] = {dPhi1,dPhi2,ntracksCharged};
1821 fMapThree[step]->Fill(propThree,strangeWeightEventAxis*strangeWeightOthers*strangeWeightThird );
1822
1823
1824 }// end of three particle correlation loop
1825
1826 }// if fThreeParticleCorr is set to true
1827
1828 }// end of inner track loop
1829
1830 } //end of outer track loop
1831
1832 fTriggerNch[step]->Fill(ntracksCharged,highPtTracksInnerEta);
1833 fTriggerNchSeeds[step]->Fill(ntracksCharged,highPtTracksInnerEta);
1834 fTriggerTracklet[step]->Fill(ntracklets);
1835
1836
1837 }//if there is at least one high pt track
1838
1839
1840 if(pindexInnerEta){// clean up array memory used for TMath::Sort
1841 delete[] pindexInnerEta;
1842 pindexInnerEta=0;
b9ad6f04 1843 }
fa8e4941 1844
1845 if(ptInnerEta){// clean up array memory used for TMath::Sort
1846 delete[] ptInnerEta;
1847 ptInnerEta=0;
b9ad6f04 1848 }
fa8e4941 1849
117f99e3 1850}
1851
7099c89a 1852
1853
117f99e3 1854//________________________________________________________________________
1855void AliAnalysisTaskMinijet::Terminate(Option_t*)
1856{
fa8e4941 1857 //terminate function is called at the end
1858 //can be used to draw histograms etc.
1859
1860
117f99e3 1861}
1862
1863//________________________________________________________________________
a640cfb3 1864Bool_t AliAnalysisTaskMinijet::SelectParticlePlusCharged(const Short_t charge, const Int_t pdg, Bool_t prim)
b9ad6f04 1865{
fa8e4941 1866 //selection of mc particle
1867 //fSelectParticles=0: use charged primaries and pi0 and k0
1868 //fSelectParticles=1: use only charged primaries
1869 //fSelectParticles=2: use only pi0 and k0
1870
1871 if(fSelectParticles==0 || fSelectParticles==2){ // in case of 2: need to count also charged particles
1872 if(charge==0){
1873 if(!(pdg==111||pdg==130||pdg==310))
1874 return false;
1875 }
1876 else{// charge !=0
1877 if(!prim)
1878 return false;
1879 }
b9ad6f04 1880 }
fa8e4941 1881
1882 else if(fSelectParticles==1){
1883 if (charge==0 || !prim){
1884 return false;
1885 }
b9ad6f04 1886 }
fa8e4941 1887
1888 else{
28cc9dd5 1889 AliInfo("Error: wrong selection of charged/pi0/k0");
fa8e4941 1890 return 0;
b9ad6f04 1891 }
fa8e4941 1892
1893 return true;
1894
b9ad6f04 1895}
1896
1897//________________________________________________________________________
a640cfb3 1898Bool_t AliAnalysisTaskMinijet::SelectParticle(const Short_t charge, const Int_t pdg, const Bool_t prim)
b9ad6f04 1899{
fa8e4941 1900 //selection of mc particle
1901 //fSelectParticles=0: use charged primaries and pi0 and k0
1902 //fSelectParticles=1: use only charged primaries
1903 //fSelectParticles=2: use only pi0 and k0
1904
1905 if(fSelectParticles==0){
1906 if(charge==0){
1907 if(!(pdg==111||pdg==130||pdg==310))
1908 return false;
1909 }
1910 else{// charge !=0
1911 if(!prim)
1912 return false;
1913 }
b9ad6f04 1914 }
fa8e4941 1915 else if (fSelectParticles==1){
1916 if (charge==0 || !prim){
1917 return false;
1918 }
b9ad6f04 1919 }
fa8e4941 1920 else if(fSelectParticles==2){
1921 if(!(pdg==111||pdg==130||pdg==310))
1922 return false;
b9ad6f04 1923 }
fa8e4941 1924
1925 return true;
1926
b9ad6f04 1927}
1928
1929//________________________________________________________________________
a640cfb3 1930Bool_t AliAnalysisTaskMinijet::CheckEvent(const Bool_t recVertex)
117f99e3 1931{
fa8e4941 1932 // This function tests the quality of an event (ESD/AOD) (rec and/or mc part)
1933 // recVertex=false: check if Mc events and stack is there, Nmc>0
1934 // recVertex=false: " + check if there is a good, reconstructed SPD vertex
1935 // defined by |z|<fVertexCut(10cm), Contributer>0, no PileUpFromSPD(3,0,8)
7099c89a 1936
fa8e4941 1937 if(fMode==0){//esd
1938
1939 //mc
1940 if(fUseMC){
1941
1942 //mc event
1943 AliMCEvent *mcEvente = (AliMCEvent*) MCEvent();
1944 if (!mcEvente) {
1945 Error("CheckEvent", "Could not retrieve MC event");
1946 return false;
1947 }
1948
1949 //stack
1950 AliStack* stackg = MCEvent()->Stack();
1951 if(!stackg) return false;
1952 Int_t ntracksg = mcEvente->GetNumberOfTracks();
1953 if(ntracksg<0) return false;
1954
1955 //vertex
1956 AliGenEventHeader* headerg = MCEvent()->GenEventHeader();
1957 TArrayF mcVg;
1958 headerg->PrimaryVertex(mcVg);
1959 //if(TMath::Abs(mcVg[0])<1e-8 && TMath::Abs(mcVg[0])<1e-8 &&
1960 // TMath::Abs(mcVg[0])<1e-8) return false;
1961 Float_t vzMCg = mcVg[2];
1962 if(TMath::Abs(vzMCg)>fVertexZCut) return false;
1963 //hasVtxMc=true;
1964 }
1965
1966 //rec
1967 if(recVertex==true){
1968 if(fESDEvent->IsPileupFromSPD(3,0,8))return false;
1969
1970 //rec vertex
1971 const AliESDVertex* vertexESDg = fESDEvent->GetPrimaryVertex(); // uses track or SPD vertexer
1972 if(!vertexESDg) return false;
1973 fVertexCheck->Fill(vertexESDg->GetNContributors());
1974 if(vertexESDg->GetNContributors()<=0)return false;
1975 Float_t fVzg= vertexESDg->GetZ();
1976 if(TMath::Abs(fVzg)>fVertexZCut) return false;
1977
1978 //rec spd vertex
1979 const AliESDVertex *vtxSPD = fESDEvent->GetPrimaryVertexSPD();
1980 if(!vtxSPD) return false;
1981 if(vtxSPD->GetNContributors()<=0)return false;
1982 Float_t fVzSPD= vtxSPD->GetZ();
1983 if(TMath::Abs(fVzSPD)>fVertexZCut) return false;
1984
1985 }
1986 return true;
7099c89a 1987 }
fa8e4941 1988
1989
1990 else if(fMode==1){ //aod
1991
1992 if(fUseMC){
9cccc98a 1993
fa8e4941 1994 //retreive MC particles from event
1995 TClonesArray *mcArray = (TClonesArray*)fAODEvent->
1996 FindListObject(AliAODMCParticle::StdBranchName());
1997 if(!mcArray){
28cc9dd5 1998 AliInfo("No MC particle branch found");
fa8e4941 1999 return false;
2000 }
2001
2002 //mc
2003 AliAODMCHeader *aodMCheader = (AliAODMCHeader *) fAODEvent->
2004 FindListObject(AliAODMCHeader::StdBranchName());
2005 Float_t vzMC = aodMCheader->GetVtxZ();
2006 if(TMath::Abs(vzMC)>fVertexZCut) return false;
2007
2008 //hasVtxMc=true;
2009 }
2010
2011 //rec
2012 if(recVertex==true){
2013 //if(fAODEvent->IsPileupFromSPD(3,0.8)) return false;
2014
2015 AliAODVertex* vertex= (AliAODVertex*)fAODEvent->GetPrimaryVertex();
2016 if(!vertex) return false;
2017 TString vtxTitle(vertex->GetTitle());// only allow vertex from tracks, no vertexer z
2018 // Printf("vtxTitle: %s",vtxTitle.Data());
2019 //if (!(vtxTitle.Contains("VertexerTracksWithConstraint"))) return false;
2020 fVertexCheck->Fill(vertex->GetNContributors());
2021 if(vertex->GetNContributors()<=0) return false;
2022 Double_t vzAOD=vertex->GetZ();
2023 // if(TMath::Abs(vzAOD)<1e-9) return false;
2024 if(TMath::Abs(vzAOD)>fVertexZCut) return false;
2025
2026 AliAODVertex* vertexSPD= (AliAODVertex*)fAODEvent->GetPrimaryVertexSPD();
2027 if(!vertexSPD) return false;
2028 if(vertexSPD->GetNContributors()<=0) return false;
2029 Double_t vzSPD=vertexSPD->GetZ();
2030 //if(TMath::Abs(vzSPD)<1e-9) return false;
2031 if(TMath::Abs(vzSPD)>fVertexZCut) return false;
2032
2033 //check TPC reconstruction: check for corrupted chunks
2034 //better: check TPCvertex, but this is not available yet in AODs
2035 Int_t nAcceptedTracksTPC=0;
2036 Int_t nAcceptedTracksITSTPC=0;
2037 for (Int_t iTracks = 0; iTracks < fAODEvent->GetNumberOfTracks(); iTracks++) {
2038 AliAODTrack *track = (AliAODTrack *)fAODEvent->GetTrack(iTracks);
2039 if (!track) continue;
2040 if(track->TestFilterBit(128) && TMath::Abs(track->Eta())<fEtaCut &&
2041 track->Pt()>fPtMin && track->Pt()<fPtMax)
2042 nAcceptedTracksTPC++;
2043 if(track->TestFilterBit(fFilterBit) && TMath::Abs(track->Eta())<fEtaCut &&
2044 track->Pt()>fPtMin && track->Pt()<fPtMax)
2045 nAcceptedTracksITSTPC++;
2046 }
2047 fCorruptedChunks->Fill(nAcceptedTracksTPC,nAcceptedTracksITSTPC);
2048 if(fRejectChunks){
2049 if(nAcceptedTracksTPC>fNTPC && nAcceptedTracksITSTPC==0)
2050 return false;//most likely corrupted chunk. No ITS tracks are reconstructed
2051 }
2052 fCorruptedChunksAfter->Fill(nAcceptedTracksTPC,nAcceptedTracksITSTPC);
2053
2054 //control histograms=================
2055 //tracklet loop
2056 Int_t ntrackletsAccept=0;
2057 AliAODTracklets * mult= (AliAODTracklets*)fAODEvent->GetTracklets();
2058 for(Int_t i=0;i<mult->GetNumberOfTracklets();++i){
2059 if(TMath::Abs(mult->GetDeltaPhi(i))<0.05 &&
2060 TMath::Abs(TMath::Log(TMath::Tan(0.5 * mult->GetTheta(i))))<fEtaCut) ++ntrackletsAccept;
2061 }
2062 Int_t nAcceptedTracks=0;
2063 for (Int_t iTracks = 0; iTracks < fAODEvent->GetNumberOfTracks(); iTracks++) {
2064 AliAODTrack *track = (AliAODTrack *)fAODEvent->GetTrack(iTracks);
2065 if (!track) continue;
2066 if(track->TestFilterBit(fFilterBit) && TMath::Abs(track->Eta())<fEtaCut
2067 && track->Pt()>fPtMin && track->Pt()<fPtMax) nAcceptedTracks++;
2068 }
2069 fNContrNtracklets->Fill(ntrackletsAccept,vertexSPD->GetNContributors());
2070 fNContrNtracks->Fill(nAcceptedTracks,vertexSPD->GetNContributors());
2071 //====================================
2072 }
2073 return true;
2074
7099c89a 2075 }
fa8e4941 2076
2077 else {
2078 Printf("Wrong mode!");
2079 return false;
7099c89a 2080 }
fa8e4941 2081
7099c89a 2082}
2083
2084//_____________________________________________________________________________
fa8e4941 2085const Double_t * AliAnalysisTaskMinijet::CreateLogAxis(const Int_t nbins,
2086 const Double_t xmin,
2087 const Double_t xmax)
7099c89a 2088{
fa8e4941 2089 // returns pointer to an array with can be used to build a logarithmic axis
2090 // it is user responsibility to delete the array
2091
2092 Double_t logxmin = TMath::Log10(xmin);
2093 Double_t logxmax = TMath::Log10(xmax);
2094 Double_t binwidth = (logxmax-logxmin)/nbins;
2095
2096 Double_t *xbins = new Double_t[nbins+1];
2097
2098 xbins[0] = xmin;
2099 for (Int_t i=1;i<=nbins;i++) {
2100 xbins[i] = xmin + TMath::Power(10,logxmin+i*binwidth);
2101 }
2102
2103 return xbins;
117f99e3 2104}
b9ad6f04 2105
7099c89a 2106//_____________________________________________________________________________
a640cfb3 2107Bool_t AliAnalysisTaskMinijet::CheckLikeSign(const Short_t chargeEventAxis,
fa8e4941 2108 const Short_t chargeOthers)
7099c89a 2109{
fa8e4941 2110 // compute if charge of two particles/tracks has same sign or different sign
2111
2112 if(fDebug>2) Printf("Charge1=%d, Charge2=%d",chargeEventAxis,chargeOthers);
2113
2114 if(chargeEventAxis<0){
2115 if(chargeOthers<0){
2116 return true;
2117 }
2118 else if(chargeOthers>0){
2119 return false;
2120 }
7099c89a 2121 }
fa8e4941 2122
2123 else if(chargeEventAxis>0){
2124 if(chargeOthers>0){
2125 return true;
2126 }
2127 else if(chargeOthers<0){
2128 return false;
2129 }
7099c89a 2130 }
fa8e4941 2131
2132 else{
28cc9dd5 2133 AliInfo("Error: Charge not lower nor higher as zero");
fa8e4941 2134 return false;
7099c89a 2135 }
fa8e4941 2136
28cc9dd5 2137 AliInfo("Error: Check values of Charge");
7099c89a 2138 return false;
7099c89a 2139}
fa8e4941 2140
2141