]> 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
variable type change (Emilia)
[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 {
179 nbinsCentr = 100;
180 minbinCentr=0;
181 maxbinCentr=100;
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)
681 ntracks = TMath::Nint(centrality);
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];
1239 Bool_t isDca= track->PropagateToDCA(fAODEvent->GetPrimaryVertex(),fAODEvent->GetMagneticField(),9999.,d0rphiz,covd0);
1240 fPropagateDca->Fill(Int_t(isDca));
1241 if(TMath::Abs(save - track->Pt())>1e-6) Printf("Before pt=%f, After pt=%f",save, track->Pt());
1242
1243 if(track->TestFilterBit(fFilterBit) && TMath::Abs(track->Eta())<fEtaCut
1244 && track->Pt()>fPtMin && track->Pt()<fPtMax){
1245
1246 nAcceptedTracks++;
1247
1248 // Printf("dca= %f", track->DCA());
1249 //save track properties in vector
1250 ptArray.push_back(track->Pt());
1251 etaArray.push_back(track->Eta());
1252 phiArray.push_back(track->Phi());
1253 chargeArray.push_back(track->Charge());
1254 fHistPt->Fill(track->Pt());
1255
1256
1257 //correction for underestimation of strangeness in Monte Carlos -> underestimation of contamination
1258 Float_t factor=1.;
1259 if(fUseMC && fCorrStrangeness && step==7){
1260 if(vtrack->GetLabel()>0){
1261 AliAODMCParticle* mcprong =(AliAODMCParticle*)mcArray->At(vtrack->GetLabel());
1262 if(mcprong){
1263 Int_t labmom = mcprong->GetMother();
1264 if(labmom>=0){
1265 AliAODMCParticle* mcmother=(AliAODMCParticle*)mcArray->At(labmom);
1266 Int_t pdgMother=0;
1267 if(mcmother) {
1268 pdgMother = mcmother->GetPdgCode();
1269 if(TMath::Abs(pdgMother)==310 || TMath::Abs(pdgMother)==130 || TMath::Abs(pdgMother)==321){ //K^0_S, K^0_L, K^+-
1270 if(track->Pt()<=1)factor=1./0.7; // values from strangeness publication
1271 else factor=1./0.6;// values from strangeness publication
1272 }
1273 if(TMath::Abs(pdgMother)==3122) { //Lambda
1274 factor=1./0.25; // values from strangeness publication
1275 }
1276 }
1277 }
1278 }
1279 }
1280 }
1281 nAcceptedTracksStrange+=factor;
1282 strangeArray.push_back(factor);
1283 fDcaXY[step]->Fill(d0rphiz[0], factor);
1284 fDcaZ[step]->Fill(d0rphiz[0], factor);
1285
1286 }
7099c89a 1287 }
fa8e4941 1288 //need to check this option for MC
1289 if(nAcceptedTracks==0) return 0;
7099c89a 1290
2b52b282 1291
fa8e4941 1292 //tracklet loop
1293 Int_t ntrackletsAccept=0;
1294 AliAODTracklets * mult= (AliAODTracklets*)fAODEvent->GetTracklets();
1295 for(Int_t i=0;i<mult->GetNumberOfTracklets();++i){
1296 if(TMath::Abs(mult->GetDeltaPhi(i))<0.05 && TMath::Abs(TMath::Log(TMath::Tan(0.5 * mult->GetTheta(i))))<fEtaCut){
1297 ++ntrackletsAccept;
1298 }
7099c89a 1299 }
fa8e4941 1300
1301
1302 nTracksTracklets.push_back(nAcceptedTracks);
1303 nTracksTracklets.push_back(ntrackletsAccept);
1304 nTracksTracklets.push_back(nAcceptedTracks);//in order to be compatible with mc analysis
1305 //where here also neutral particles are counted.
1306
1307
1308 fVzEvent= vzAOD;
1309 if(step==5 || step==2){
1310 fNRecAccept = nAcceptedTracks; // needed for MC case //step5 = TrigVtxRecNrec
1311 fNRecAcceptTracklet = ntrackletsAccept; // needed for MC case //step5 = TrigVtxRecNrec
7099c89a 1312 }
fa8e4941 1313 if(step==7)fNRecAcceptStrangeCorr = nAcceptedTracksStrange;
1314
1315 return fNRecAccept; // at the moment, always return reconstructed multiplicity
1316
7099c89a 1317}
1318
1319//________________________________________________________________________
11722f75 1320Double_t AliAnalysisTaskMinijet::ReadEventAODRecMcProp( vector<Float_t> &ptArray, vector<Float_t> &etaArray,
fa8e4941 1321 vector<Float_t> &phiArray, vector<Short_t> &chargeArray,
1322 vector<Float_t> &strangeArray,
11722f75 1323 vector<Double_t> &nTracksTracklets, const Int_t step)
7099c89a 1324{
fa8e4941 1325 // gives back the number of AOD tracks and pointer to arrays with track
1326 // properties (pt, eta, phi)
1327
1328
1329 ptArray.clear();
1330 etaArray.clear();
1331 phiArray.clear();
1332 chargeArray.clear();
1333 strangeArray.clear();
1334 nTracksTracklets.clear();
1335
1336
1337 // Retreive the number of tracks for this event.
1338 Int_t ntracks = fAODEvent->GetNumberOfTracks();
1339 if(fDebug>1) Printf("AOD tracks: %d", ntracks);
1340
1341
1342 //get array of mc particles
1343 TClonesArray *mcArray = (TClonesArray*)fAODEvent->
7099c89a 1344 FindListObject(AliAODMCParticle::StdBranchName());
fa8e4941 1345 if(!mcArray){
28cc9dd5 1346 AliInfo("No MC particle branch found");
fa8e4941 1347 return kFALSE;
117f99e3 1348 }
fa8e4941 1349
1350 AliAODVertex* vtx= (AliAODVertex*)fAODEvent->GetPrimaryVertexSPD();//GetPrimaryVertex()
1351 Double_t vzAOD=vtx->GetZ();
1352 fVertexZ[step]->Fill(vzAOD);
1353
11722f75 1354 Double_t nAcceptedTracks=0;
fa8e4941 1355 for (Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
1356 AliAODTrack *track = (AliAODTrack *)fAODEvent->GetTrack(iTracks);
1357
1358 AliVParticle *vtrack = fAODEvent->GetTrack(iTracks);
1359
1360 if (!track) {
1361 Error("ReadEventAODRecMcProp", "Could not receive track %d", iTracks);
1362 continue;
1363 }
1364
1365 //use only tracks from primaries
1366 if(fAnalysePrimOnly){
1367 if(vtrack->GetLabel()<=0)continue;
1368 if(!(static_cast<AliAODMCParticle*>(mcArray->At(vtrack->GetLabel()))->IsPhysicalPrimary()))continue;
1369 }
1370
1371 if(track->TestFilterBit(fFilterBit) && TMath::Abs(track->Eta())<fEtaCut &&
1372 track->Pt()>fPtMin && track->Pt()<fPtMax){
1373
1374 nAcceptedTracks++;
1375 Float_t factor=1.;
1376
1377 //save track properties in vector
1378 if(vtrack->GetLabel()<=0){ //fake tracks before "label<0", but crash in AOD079 // what is the meaning of label 0
1379 // Printf("Fake track");
1380 // continue;
1381 ptArray.push_back(track->Pt());
1382 etaArray.push_back(track->Eta());
1383 phiArray.push_back(track->Phi());
1384 chargeArray.push_back(track->Charge());
2b52b282 1385
fa8e4941 1386 }
1387 else{//mc properties
1388 AliAODMCParticle *partOfTrack = (AliAODMCParticle*)mcArray->At(vtrack->GetLabel());
1389 if(!partOfTrack) Printf("label=%d", vtrack->GetLabel());
1390 if(partOfTrack){
1391 ptArray.push_back(partOfTrack->Pt());
1392 etaArray.push_back(partOfTrack->Eta());
1393 phiArray.push_back(partOfTrack->Phi());
1394 chargeArray.push_back(vtrack->Charge());//partOfTrack?
1395
1396 //correction for underestimation of strangeness in Monte Carlos -> underestimation of contamination
1397 if(fUseMC && fCorrStrangeness && step==6){
1398 Int_t labmom = partOfTrack->GetMother();
1399 if(labmom>=0){
1400 AliAODMCParticle* mcmother=(AliAODMCParticle*)mcArray->At(labmom);
1401 Int_t pdgMother=0;
1402 if(mcmother) {
1403 pdgMother = mcmother->GetPdgCode();
1404 if(TMath::Abs(pdgMother)==310 || TMath::Abs(pdgMother)==130 || TMath::Abs(pdgMother)==321){ //K^0_S, K^0_L, K^+-
1405 if(track->Pt()<=1)factor=1./0.7; // values from strangeness publication
1406 else factor=1./0.6;// values from strangeness publication
1407 }
1408 if(TMath::Abs(pdgMother)==3122) { //Lambda
1409 factor=1./0.25; // values from strangeness publication
1410 }
1411 }
1412 }
1413 }
1414 }
1415 }
1416 strangeArray.push_back(factor);
1417
1418 }
7099c89a 1419 }
fa8e4941 1420 //need to check this option for MC
1421 if(nAcceptedTracks==0) return 0;
1422
1423 //tracklet loop
1424 Int_t ntrackletsAccept=0;
1425 AliAODTracklets * mult= (AliAODTracklets*)fAODEvent->GetTracklets();
1426 for(Int_t i=0;i<mult->GetNumberOfTracklets();++i){
1427 if(TMath::Abs(mult->GetDeltaPhi(i))<0.05 && TMath::Abs(TMath::Log(TMath::Tan(0.5 * mult->GetTheta(i))))<fEtaCut ){
1428 ++ntrackletsAccept;
1429 }
117f99e3 1430 }
fa8e4941 1431
1432
1433 nTracksTracklets.push_back(nAcceptedTracks);
1434 nTracksTracklets.push_back(ntrackletsAccept);
1435 nTracksTracklets.push_back(nAcceptedTracks);//in order to be compatible with mc analysis
1436 //where here also neutral particles are counted.
1437
1438
1439 //check vertex (mc)
1440 AliAODMCHeader *aodMCheader = (AliAODMCHeader *) fAODEvent->
7099c89a 1441 FindListObject(AliAODMCHeader::StdBranchName());
fa8e4941 1442 Float_t vzMC = aodMCheader->GetVtxZ();
1443
1444 fVzEvent= vzMC;
1445 return fNRecAccept;//this is the rec value from step 5
1446
1447}
7099c89a 1448
1449
117f99e3 1450
1451//________________________________________________________________________
11722f75 1452Double_t AliAnalysisTaskMinijet::ReadEventAODMC( vector<Float_t> &ptArray, vector<Float_t> &etaArray,
fa8e4941 1453 vector<Float_t> &phiArray, vector<Short_t> &chargeArray,
1454 vector<Float_t> &strangeArray,
11722f75 1455 vector<Double_t> &nTracksTracklets, const Int_t step)
117f99e3 1456{
fa8e4941 1457 // gives back the number of AOD MC particles and pointer to arrays with particle
1458 // properties (pt, eta, phi)
1459
1460 ptArray.clear();
1461 etaArray.clear();
1462 phiArray.clear();
1463 chargeArray.clear();
1464 strangeArray.clear();
1465 nTracksTracklets.clear();
1466
1467 //check vertex
1468 AliAODMCHeader *aodMCheader = (AliAODMCHeader *) fAODEvent->
7099c89a 1469 FindListObject(AliAODMCHeader::StdBranchName());
fa8e4941 1470 Float_t vzMC = aodMCheader->GetVtxZ();
1471 if(step==1){
1472 fVertexZ[0]->Fill(vzMC);
1473 }
1474 fVertexZ[step]->Fill(vzMC);
1475
1476 //retreive MC particles from event
1477 TClonesArray *mcArray = (TClonesArray*)fAODEvent->
b7c0438e 1478 FindListObject(AliAODMCParticle::StdBranchName());
fa8e4941 1479 if(!mcArray){
28cc9dd5 1480 AliInfo("No MC particle branch found");
fa8e4941 1481 return kFALSE;
117f99e3 1482 }
fa8e4941 1483
1484 Int_t ntracks = mcArray->GetEntriesFast();
1485 if(fDebug>1)Printf("MC particles: %d", ntracks);
1486
1487
1488 // Track loop: chek how many particles will be accepted
1489 //Float_t vzMC=0.;
11722f75 1490 Double_t nChargedPrim=0;
1491 Double_t nAllPrim=0;
1492 Double_t nPseudoTracklets=0;
fa8e4941 1493 for (Int_t it = 0; it < ntracks; it++) {
1494 AliAODMCParticle *track = (AliAODMCParticle*)mcArray->At(it);
1495 if (!track) {
1496 Error("ReadEventAODMC", "Could not receive particle %d", it);
1497 continue;
1498 }
1499
1500 if(!SelectParticlePlusCharged(
1501 track->Charge(),
1502 track->PdgCode(),
1503 track->IsPhysicalPrimary()
1504 )
1505 )
1506 continue;
1507
1508 if(TMath::Abs(track->Eta())<fEtaCut && track->Pt()>0.0)++nPseudoTracklets; //0.035
1509 if(TMath::Abs(track->Eta())>fEtaCut || track->Pt()<fPtMin || track->Pt()>fPtMax) continue;
1510
1511 nAllPrim++;
1512 if(track->Charge()!=0) nChargedPrim++;
1513
b7c0438e 1514 }
fa8e4941 1515
1516
1517 if(nAllPrim==0) return 0;
1518 if(nChargedPrim==0) return 0;
1519
1520 //generate array with size of number of accepted tracks
1521 fChargedPi0->Fill(nAllPrim,nChargedPrim);
1522
1523
1524 // Track loop: fill arrays for accepted tracks
1525 // Int_t iChargedPrim=0;
1526 for (Int_t it = 0; it < ntracks; it++) {
1527 AliAODMCParticle *track = (AliAODMCParticle*)mcArray->At(it);
1528 if (!track) {
1529 Error("ReadEventAODMC", "Could not receive particle %d", it);
1530 continue;
1531 }
1532
1533 if(!SelectParticle(
1534 track->Charge(),
1535 track->PdgCode(),
1536 track->IsPhysicalPrimary()
1537 )
1538 )
1539 continue;
1540
1541 if(TMath::Abs(track->Eta())>fEtaCut || track->Pt()<fPtMin || track->Pt()>fPtMax) continue;
1542
1543 fHistPtMC->Fill(track->Pt());
1544 ptArray.push_back(track->Pt());
1545 etaArray.push_back(track->Eta());
1546 phiArray.push_back(track->Phi());
1547 chargeArray.push_back(track->Charge());
1548 strangeArray.push_back(1);
1549 }
1550
1551 nTracksTracklets.push_back(nChargedPrim);
1552 nTracksTracklets.push_back(nPseudoTracklets);
1553 if(fSelectParticles!=2){
1554 nTracksTracklets.push_back(nAllPrim);
1555 }
1556 else{
1557 nTracksTracklets.push_back(nAllPrim-nChargedPrim); // only neutral
1558 }
1559
1560
1561
1562 fVzEvent= vzMC;
1563 fNMcPrimAccept = nChargedPrim;
1564 fNMcPrimAcceptTracklet = nPseudoTracklets;
1565
1566 if(step==1){ // step 1 = Trig All Mc Nmc
1567 fNmcNch->Fill( fNMcPrimAccept,fNRecAccept);
1568 fPNmcNch->Fill(fNMcPrimAccept,fNRecAccept);
1569 fNmcNchTracklet->Fill( fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
1570 fPNmcNchTracklet->Fill(fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
1571 }
1572 if(step==3){ // step 3 = Trig vtx Mc
1573 fNmcNchVtx->Fill( fNMcPrimAccept,fNRecAccept);
1574 fNmcNchVtxStrangeCorr->Fill( fNMcPrimAccept,fNRecAcceptStrangeCorr);
1575 fPNmcNchVtx->Fill(fNMcPrimAccept,fNRecAccept);
1576 fNmcNchVtxTracklet->Fill( fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
1577 fPNmcNchVtxTracklet->Fill(fNMcPrimAcceptTracklet,fNRecAcceptTracklet);
1578 }
1579 return fNRecAccept; // rec value from step 5 or step 2
1580
1581
117f99e3 1582}
1583
1584//________________________________________________________________________
fa8e4941 1585void AliAnalysisTaskMinijet::Analyse(const vector<Float_t> &pt,
1586 const vector<Float_t> &eta,
1587 const vector<Float_t> &phi,
1588 const vector<Short_t> &charge,
1589 const vector<Float_t> &strangeWeight,
11722f75 1590 const Double_t ntracksCharged,
fa8e4941 1591 const Int_t ntracklets,
1592 const Int_t nAll,
1593 const Int_t step)
117f99e3 1594{
fa8e4941 1595
1596 // analyse track properties (comming from either ESDs or AODs) in order to compute
1597 // mini jet activity (chared tracks) as function of charged multiplicity
1598
1599 fStep->Fill(step);
1600
1601 if(fDebug){
1602 Printf("Analysis Step=%d", step);
1603 if(fDebug>2){
1604 Printf("nAll=%d",nAll);
11722f75 1605 Printf("nCharged=%f",ntracksCharged);
fa8e4941 1606 for (Int_t i = 0; i < nAll; i++) {
1607 Printf("pt[%d]=%f",i,pt[i]);
1608 }
1609 }
7099c89a 1610 }
fa8e4941 1611
1612 //calculation of mean pt for all tracks in case of step==0
1613 if(step==5 || step ==2){ // rec step
1614 Double_t meanPt=0.;
1615 Double_t leadingPt=0.;
1616 for (Int_t i = 0; i < nAll; i++) {
1617 if(pt[i]<0.01)continue;
1618 meanPt+=pt[i];
1619 if(leadingPt<pt[i])leadingPt=pt[i];
1620 }
1621 meanPt=meanPt/nAll;
1622 fMeanPtRec=meanPt;
1623 fLeadingPtRec=leadingPt;
1624 }
1625
1626 Double_t propEvent[4] = {ntracksCharged,fVzEvent,fMeanPtRec, fLeadingPtRec}; //vz: {rec, mc, mc}, meanPt and Nrec is always rec value
1627 fMapEvent[step]->Fill(propEvent);
1628
1629 fNcharge[step]->Fill(ntracksCharged);
1630
1631 Float_t ptEventAxis=0; // pt event axis
1632 Float_t etaEventAxis=0; // eta event axis
1633 Float_t phiEventAxis=0; // phi event axis
1634 Short_t chargeEventAxis=0; // charge event axis
1635 Float_t strangeWeightEventAxis=0; // strange weight event axis
1636
1637 Float_t ptOthers = 0; // pt others // for all other tracks around event axis -> see loop
1638 Float_t etaOthers = 0; // eta others
1639 Float_t phiOthers = 0; // phi others
1640 Short_t chargeOthers = 0; // charge others
1641 Float_t strangeWeightOthers = 0; // strange weight others
1642
1643 Int_t *pindexInnerEta = new Int_t [nAll+1];
1644 Float_t *ptInnerEta = new Float_t[nAll+1];
1645
1646
1647
7099c89a 1648 for (Int_t i = 0; i < nAll; i++) {
fa8e4941 1649
1650 if(pt[i]<0.01)continue;
1651
1652 //fill single particle correction for first step of pair correction
1653 Double_t propAll[3] = {eta[i],pt[i],ntracksCharged };
1654 fMapAll[step]->Fill(propAll, strangeWeight[i]);
1655
1656
1657 //filling of simple check plots
1658 if(pt[i]<0.7) continue;
1659 fPt[step] -> Fill( pt[i]);
1660 fEta[step] -> Fill(eta[i]);
1661 fPhi[step] -> Fill(phi[i]);
1662 fPhiEta[step]-> Fill(phi[i], eta[i]);
1663
1664 pindexInnerEta[i]=0; //set all values to zero
1665 //fill new array for eta check
1666 ptInnerEta[i]= pt[i];
1667
7099c89a 1668 }
fa8e4941 1669
1670
1671
1672 // define event axis: leading or random track (with pt>fTriggerPtCut)
1673 // ---------------------------------------
1674 Int_t highPtTracks=0;
1675 Int_t highPtTracksInnerEta=0;
1676 // Double_t highPtTracksInnerEtaCorr=0;
1677 Int_t mult09=0;
1678
1679 //count high pt tracks and high pt tracks in acceptance for seeds
1680 for(Int_t i=0;i<nAll;i++){
1681
1682 if(pt[i]<0.01)continue;
1683
1684 if(TMath::Abs(eta[i])<0.9){
1685 mult09++;
1686 }
1687
1688 if(pt[i]>fTriggerPtCut) {
1689 highPtTracks++;
1690 }
1691
1692 // seed should be place in middle of acceptance, that complete cone is in acceptance
1693 if(pt[i]>fTriggerPtCut && TMath::Abs(eta[i])<fEtaCutSeed && charge[i]!=0){
1694
1695 // Printf("eta=%f", eta[i]);
1696 highPtTracksInnerEta++;
1697
1698 }
1699 else{
1700 ptInnerEta[i]=0;
1701 }
9ed461df 1702 }
fa8e4941 1703
1704
1705 //sort array in order to get highest pt tracks first
1706 //index can be computed with pindexInnerEta[number]
1707 if(nAll) TMath::Sort(nAll, ptInnerEta, pindexInnerEta, kTRUE);
1708
1709 // plot of multiplicity distributions
1710 fNch07Nch[step]->Fill(ntracksCharged, highPtTracksInnerEta);
1711 fPNch07Nch[step]->Fill(ntracksCharged, highPtTracksInnerEta);
1712
1713 if(ntracklets){
1714 fNch07Tracklet[step]->Fill(ntracklets, highPtTracksInnerEta);//only counts tracks which can be used as seeds
1715 fNchTracklet[step]->Fill(ntracklets, ntracksCharged);
1716 fPNch07Tracklet[step]->Fill(ntracklets, highPtTracksInnerEta);//only counts tracks which can be used as seeds
b9ad6f04 1717 }
fa8e4941 1718
1719 //analysis can only be performed with event axis, defined by high pt track
1720
1721
1722 if(highPtTracks>0 && highPtTracksInnerEta>0){
1723
1724 // build pairs in two track loops
1725 // loop over all possible trigger particles (defined by pt_trig and eta_acceptance)
1726 for(Int_t axis=0;(axis<nAll) && (pt[pindexInnerEta[axis]]>=fTriggerPtCut); axis++){
1727
1728 //EventAxisRandom track properties
1729 ptEventAxis = pt [pindexInnerEta[axis]];
1730 etaEventAxis = eta[pindexInnerEta[axis]];
1731 phiEventAxis = phi[pindexInnerEta[axis]];
1732 chargeEventAxis = charge[pindexInnerEta[axis]];
1733 strangeWeightEventAxis = strangeWeight[pindexInnerEta[axis]];
1734 fPtSeed[step] -> Fill( ptEventAxis);
1735 fEtaSeed[step] -> Fill(etaEventAxis);
1736 fPhiSeed[step] -> Fill(phiEventAxis);
1737
1738
1739 Double_t prop[3] = {etaEventAxis,ptEventAxis,ntracksCharged };
1740 fMapSingleTrig[step]->Fill(prop, strangeWeightEventAxis);
1741
1742 //associated tracks
1743 for (Int_t iTrack = axis+1; iTrack < nAll; iTrack++) {
1744
1745 if(pt[pindexInnerEta[iTrack]]<fAssociatePtCut) continue;
1746
1747 if(fSelectParticlesAssoc==1){
1748 if(charge[pindexInnerEta[iTrack]]==0)continue;
1749 }
1750 if(fSelectParticlesAssoc==2){
1751 if(charge[pindexInnerEta[iTrack]]!=0)continue;
1752 }
1753
1754
1755 ptOthers = pt [pindexInnerEta[iTrack]];
1756 etaOthers = eta[pindexInnerEta[iTrack]];
1757 phiOthers = phi[pindexInnerEta[iTrack]];
1758 chargeOthers = charge[pindexInnerEta[iTrack]];
1759 strangeWeightOthers = strangeWeight[pindexInnerEta[iTrack]];
1760
1761
1762 //plot only properties of tracks with pt>ptassoc
1763 fPtOthers[step] -> Fill( ptOthers);
1764 fEtaOthers[step] -> Fill(etaOthers);
1765 fPhiOthers[step] -> Fill(phiOthers);
1766 fPtEtaOthers[step] -> Fill(ptOthers, etaOthers);
1767
1768 // if(fDebug>2)Printf("%f, %f", pt[pindexInnerEta[axis]], pt[pindexInnerEta[iTrack]]);
1769
1770 Float_t dPhi = phiOthers-phiEventAxis;
1771 if(dPhi>1.5*TMath::Pi()) dPhi = dPhi-2*TMath::Pi();
1772 else if(dPhi<-0.5*TMath::Pi())dPhi=dPhi+2*TMath::Pi();
1773 Float_t dEta=etaOthers-etaEventAxis;
1774
1775
1776 fDPhiDEtaEventAxis[step]->Fill(dPhi, dEta, strangeWeightEventAxis*strangeWeightOthers);
1777 fDPhiEventAxis[step]->Fill(dPhi, strangeWeightEventAxis*strangeWeightOthers);
1778
1779 //check outliers
28cc9dd5 1780 if(ptEventAxis< 0.4 || ptEventAxis > 100) AliInfo("particles out of range pt");
1781 if(ntracksCharged<-1 || ntracksCharged>1500) AliInfo("particles out of range ncharge");
fa8e4941 1782 if(TMath::Abs(dEta)>2*fEtaCut) {
28cc9dd5 1783 AliInfo("particles out of range dEta");
c8ade419 1784 AliInfo(Form("eta1=%f, eta2=%f", etaOthers, etaEventAxis));
1785 AliInfo(Form("step=%d",step));
fa8e4941 1786 }
1787 if(dPhi<-0.5*TMath::Pi() || dPhi>1.5*TMath::Pi()){
c8ade419 1788 AliInfo(Form("particles out of range dPhi"));
1789 AliInfo(Form("phi1=%f, phi2=%f", phiOthers, phiEventAxis));
fa8e4941 1790 }
1791
1792 Bool_t isLikeSign = CheckLikeSign(chargeEventAxis, chargeOthers);
1793
1794 Double_t prop6[6] = {ptEventAxis,ptOthers,dEta,dPhi,ntracksCharged, isLikeSign };
1795 fMapPair[step]->Fill(prop6, strangeWeightEventAxis*strangeWeightOthers);
1796
1797 //thrid track loop (Andreas: three particle auto-correlations)
1798 if(fThreeParticleCorr){
1799 for (Int_t third = iTrack+1; third < nAll; third++) {
1800 if(pt[pindexInnerEta[iTrack]]<fTriggerPtCut) continue;
1801 if(pt[pindexInnerEta[third]]<fTriggerPtCut) continue;
1802
1803 //dphi1
1804 Float_t dPhi1 = phiEventAxis - phiOthers;
1805 if(dPhi1>1.5*TMath::Pi()) dPhi1 = dPhi1-2*TMath::Pi();
1806 else if(dPhi1<-0.5*TMath::Pi())dPhi1=dPhi1+2*TMath::Pi();
1807
1808 Float_t phiThird = phi[pindexInnerEta[third]];
1809 Float_t strangeWeightThird = strangeWeight[pindexInnerEta[third]];
1810
1811 //dphi2
1812 Float_t dPhi2 = phiEventAxis - phiThird;
1813 if(dPhi2>1.5*TMath::Pi()) dPhi2 = dPhi2-2*TMath::Pi();
1814 else if(dPhi2<-0.5*TMath::Pi())dPhi2=dPhi2+2*TMath::Pi();
1815
1816 fDPhi1DPhi2[step]-> Fill(dPhi1, dPhi2, strangeWeightEventAxis*strangeWeightOthers*strangeWeightThird);
1817 Double_t propThree[3] = {dPhi1,dPhi2,ntracksCharged};
1818 fMapThree[step]->Fill(propThree,strangeWeightEventAxis*strangeWeightOthers*strangeWeightThird );
1819
1820
1821 }// end of three particle correlation loop
1822
1823 }// if fThreeParticleCorr is set to true
1824
1825 }// end of inner track loop
1826
1827 } //end of outer track loop
1828
1829 fTriggerNch[step]->Fill(ntracksCharged,highPtTracksInnerEta);
1830 fTriggerNchSeeds[step]->Fill(ntracksCharged,highPtTracksInnerEta);
1831 fTriggerTracklet[step]->Fill(ntracklets);
1832
1833
1834 }//if there is at least one high pt track
1835
1836
1837 if(pindexInnerEta){// clean up array memory used for TMath::Sort
1838 delete[] pindexInnerEta;
1839 pindexInnerEta=0;
b9ad6f04 1840 }
fa8e4941 1841
1842 if(ptInnerEta){// clean up array memory used for TMath::Sort
1843 delete[] ptInnerEta;
1844 ptInnerEta=0;
b9ad6f04 1845 }
fa8e4941 1846
117f99e3 1847}
1848
7099c89a 1849
1850
117f99e3 1851//________________________________________________________________________
1852void AliAnalysisTaskMinijet::Terminate(Option_t*)
1853{
fa8e4941 1854 //terminate function is called at the end
1855 //can be used to draw histograms etc.
1856
1857
117f99e3 1858}
1859
1860//________________________________________________________________________
a640cfb3 1861Bool_t AliAnalysisTaskMinijet::SelectParticlePlusCharged(const Short_t charge, const Int_t pdg, Bool_t prim)
b9ad6f04 1862{
fa8e4941 1863 //selection of mc particle
1864 //fSelectParticles=0: use charged primaries and pi0 and k0
1865 //fSelectParticles=1: use only charged primaries
1866 //fSelectParticles=2: use only pi0 and k0
1867
1868 if(fSelectParticles==0 || fSelectParticles==2){ // in case of 2: need to count also charged particles
1869 if(charge==0){
1870 if(!(pdg==111||pdg==130||pdg==310))
1871 return false;
1872 }
1873 else{// charge !=0
1874 if(!prim)
1875 return false;
1876 }
b9ad6f04 1877 }
fa8e4941 1878
1879 else if(fSelectParticles==1){
1880 if (charge==0 || !prim){
1881 return false;
1882 }
b9ad6f04 1883 }
fa8e4941 1884
1885 else{
28cc9dd5 1886 AliInfo("Error: wrong selection of charged/pi0/k0");
fa8e4941 1887 return 0;
b9ad6f04 1888 }
fa8e4941 1889
1890 return true;
1891
b9ad6f04 1892}
1893
1894//________________________________________________________________________
a640cfb3 1895Bool_t AliAnalysisTaskMinijet::SelectParticle(const Short_t charge, const Int_t pdg, const Bool_t prim)
b9ad6f04 1896{
fa8e4941 1897 //selection of mc particle
1898 //fSelectParticles=0: use charged primaries and pi0 and k0
1899 //fSelectParticles=1: use only charged primaries
1900 //fSelectParticles=2: use only pi0 and k0
1901
1902 if(fSelectParticles==0){
1903 if(charge==0){
1904 if(!(pdg==111||pdg==130||pdg==310))
1905 return false;
1906 }
1907 else{// charge !=0
1908 if(!prim)
1909 return false;
1910 }
b9ad6f04 1911 }
fa8e4941 1912 else if (fSelectParticles==1){
1913 if (charge==0 || !prim){
1914 return false;
1915 }
b9ad6f04 1916 }
fa8e4941 1917 else if(fSelectParticles==2){
1918 if(!(pdg==111||pdg==130||pdg==310))
1919 return false;
b9ad6f04 1920 }
fa8e4941 1921
1922 return true;
1923
b9ad6f04 1924}
1925
1926//________________________________________________________________________
a640cfb3 1927Bool_t AliAnalysisTaskMinijet::CheckEvent(const Bool_t recVertex)
117f99e3 1928{
fa8e4941 1929 // This function tests the quality of an event (ESD/AOD) (rec and/or mc part)
1930 // recVertex=false: check if Mc events and stack is there, Nmc>0
1931 // recVertex=false: " + check if there is a good, reconstructed SPD vertex
1932 // defined by |z|<fVertexCut(10cm), Contributer>0, no PileUpFromSPD(3,0,8)
7099c89a 1933
fa8e4941 1934 if(fMode==0){//esd
1935
1936 //mc
1937 if(fUseMC){
1938
1939 //mc event
1940 AliMCEvent *mcEvente = (AliMCEvent*) MCEvent();
1941 if (!mcEvente) {
1942 Error("CheckEvent", "Could not retrieve MC event");
1943 return false;
1944 }
1945
1946 //stack
1947 AliStack* stackg = MCEvent()->Stack();
1948 if(!stackg) return false;
1949 Int_t ntracksg = mcEvente->GetNumberOfTracks();
1950 if(ntracksg<0) return false;
1951
1952 //vertex
1953 AliGenEventHeader* headerg = MCEvent()->GenEventHeader();
1954 TArrayF mcVg;
1955 headerg->PrimaryVertex(mcVg);
1956 //if(TMath::Abs(mcVg[0])<1e-8 && TMath::Abs(mcVg[0])<1e-8 &&
1957 // TMath::Abs(mcVg[0])<1e-8) return false;
1958 Float_t vzMCg = mcVg[2];
1959 if(TMath::Abs(vzMCg)>fVertexZCut) return false;
1960 //hasVtxMc=true;
1961 }
1962
1963 //rec
1964 if(recVertex==true){
1965 if(fESDEvent->IsPileupFromSPD(3,0,8))return false;
1966
1967 //rec vertex
1968 const AliESDVertex* vertexESDg = fESDEvent->GetPrimaryVertex(); // uses track or SPD vertexer
1969 if(!vertexESDg) return false;
1970 fVertexCheck->Fill(vertexESDg->GetNContributors());
1971 if(vertexESDg->GetNContributors()<=0)return false;
1972 Float_t fVzg= vertexESDg->GetZ();
1973 if(TMath::Abs(fVzg)>fVertexZCut) return false;
1974
1975 //rec spd vertex
1976 const AliESDVertex *vtxSPD = fESDEvent->GetPrimaryVertexSPD();
1977 if(!vtxSPD) return false;
1978 if(vtxSPD->GetNContributors()<=0)return false;
1979 Float_t fVzSPD= vtxSPD->GetZ();
1980 if(TMath::Abs(fVzSPD)>fVertexZCut) return false;
1981
1982 }
1983 return true;
7099c89a 1984 }
fa8e4941 1985
1986
1987 else if(fMode==1){ //aod
1988
1989 if(fUseMC){
9cccc98a 1990
fa8e4941 1991 //retreive MC particles from event
1992 TClonesArray *mcArray = (TClonesArray*)fAODEvent->
1993 FindListObject(AliAODMCParticle::StdBranchName());
1994 if(!mcArray){
28cc9dd5 1995 AliInfo("No MC particle branch found");
fa8e4941 1996 return false;
1997 }
1998
1999 //mc
2000 AliAODMCHeader *aodMCheader = (AliAODMCHeader *) fAODEvent->
2001 FindListObject(AliAODMCHeader::StdBranchName());
2002 Float_t vzMC = aodMCheader->GetVtxZ();
2003 if(TMath::Abs(vzMC)>fVertexZCut) return false;
2004
2005 //hasVtxMc=true;
2006 }
2007
2008 //rec
2009 if(recVertex==true){
2010 //if(fAODEvent->IsPileupFromSPD(3,0.8)) return false;
2011
2012 AliAODVertex* vertex= (AliAODVertex*)fAODEvent->GetPrimaryVertex();
2013 if(!vertex) return false;
2014 TString vtxTitle(vertex->GetTitle());// only allow vertex from tracks, no vertexer z
2015 // Printf("vtxTitle: %s",vtxTitle.Data());
2016 //if (!(vtxTitle.Contains("VertexerTracksWithConstraint"))) return false;
2017 fVertexCheck->Fill(vertex->GetNContributors());
2018 if(vertex->GetNContributors()<=0) return false;
2019 Double_t vzAOD=vertex->GetZ();
2020 // if(TMath::Abs(vzAOD)<1e-9) return false;
2021 if(TMath::Abs(vzAOD)>fVertexZCut) return false;
2022
2023 AliAODVertex* vertexSPD= (AliAODVertex*)fAODEvent->GetPrimaryVertexSPD();
2024 if(!vertexSPD) return false;
2025 if(vertexSPD->GetNContributors()<=0) return false;
2026 Double_t vzSPD=vertexSPD->GetZ();
2027 //if(TMath::Abs(vzSPD)<1e-9) return false;
2028 if(TMath::Abs(vzSPD)>fVertexZCut) return false;
2029
2030 //check TPC reconstruction: check for corrupted chunks
2031 //better: check TPCvertex, but this is not available yet in AODs
2032 Int_t nAcceptedTracksTPC=0;
2033 Int_t nAcceptedTracksITSTPC=0;
2034 for (Int_t iTracks = 0; iTracks < fAODEvent->GetNumberOfTracks(); iTracks++) {
2035 AliAODTrack *track = (AliAODTrack *)fAODEvent->GetTrack(iTracks);
2036 if (!track) continue;
2037 if(track->TestFilterBit(128) && TMath::Abs(track->Eta())<fEtaCut &&
2038 track->Pt()>fPtMin && track->Pt()<fPtMax)
2039 nAcceptedTracksTPC++;
2040 if(track->TestFilterBit(fFilterBit) && TMath::Abs(track->Eta())<fEtaCut &&
2041 track->Pt()>fPtMin && track->Pt()<fPtMax)
2042 nAcceptedTracksITSTPC++;
2043 }
2044 fCorruptedChunks->Fill(nAcceptedTracksTPC,nAcceptedTracksITSTPC);
2045 if(fRejectChunks){
2046 if(nAcceptedTracksTPC>fNTPC && nAcceptedTracksITSTPC==0)
2047 return false;//most likely corrupted chunk. No ITS tracks are reconstructed
2048 }
2049 fCorruptedChunksAfter->Fill(nAcceptedTracksTPC,nAcceptedTracksITSTPC);
2050
2051 //control histograms=================
2052 //tracklet loop
2053 Int_t ntrackletsAccept=0;
2054 AliAODTracklets * mult= (AliAODTracklets*)fAODEvent->GetTracklets();
2055 for(Int_t i=0;i<mult->GetNumberOfTracklets();++i){
2056 if(TMath::Abs(mult->GetDeltaPhi(i))<0.05 &&
2057 TMath::Abs(TMath::Log(TMath::Tan(0.5 * mult->GetTheta(i))))<fEtaCut) ++ntrackletsAccept;
2058 }
2059 Int_t nAcceptedTracks=0;
2060 for (Int_t iTracks = 0; iTracks < fAODEvent->GetNumberOfTracks(); iTracks++) {
2061 AliAODTrack *track = (AliAODTrack *)fAODEvent->GetTrack(iTracks);
2062 if (!track) continue;
2063 if(track->TestFilterBit(fFilterBit) && TMath::Abs(track->Eta())<fEtaCut
2064 && track->Pt()>fPtMin && track->Pt()<fPtMax) nAcceptedTracks++;
2065 }
2066 fNContrNtracklets->Fill(ntrackletsAccept,vertexSPD->GetNContributors());
2067 fNContrNtracks->Fill(nAcceptedTracks,vertexSPD->GetNContributors());
2068 //====================================
2069 }
2070 return true;
2071
7099c89a 2072 }
fa8e4941 2073
2074 else {
2075 Printf("Wrong mode!");
2076 return false;
7099c89a 2077 }
fa8e4941 2078
7099c89a 2079}
2080
2081//_____________________________________________________________________________
fa8e4941 2082const Double_t * AliAnalysisTaskMinijet::CreateLogAxis(const Int_t nbins,
2083 const Double_t xmin,
2084 const Double_t xmax)
7099c89a 2085{
fa8e4941 2086 // returns pointer to an array with can be used to build a logarithmic axis
2087 // it is user responsibility to delete the array
2088
2089 Double_t logxmin = TMath::Log10(xmin);
2090 Double_t logxmax = TMath::Log10(xmax);
2091 Double_t binwidth = (logxmax-logxmin)/nbins;
2092
2093 Double_t *xbins = new Double_t[nbins+1];
2094
2095 xbins[0] = xmin;
2096 for (Int_t i=1;i<=nbins;i++) {
2097 xbins[i] = xmin + TMath::Power(10,logxmin+i*binwidth);
2098 }
2099
2100 return xbins;
117f99e3 2101}
b9ad6f04 2102
7099c89a 2103//_____________________________________________________________________________
a640cfb3 2104Bool_t AliAnalysisTaskMinijet::CheckLikeSign(const Short_t chargeEventAxis,
fa8e4941 2105 const Short_t chargeOthers)
7099c89a 2106{
fa8e4941 2107 // compute if charge of two particles/tracks has same sign or different sign
2108
2109 if(fDebug>2) Printf("Charge1=%d, Charge2=%d",chargeEventAxis,chargeOthers);
2110
2111 if(chargeEventAxis<0){
2112 if(chargeOthers<0){
2113 return true;
2114 }
2115 else if(chargeOthers>0){
2116 return false;
2117 }
7099c89a 2118 }
fa8e4941 2119
2120 else if(chargeEventAxis>0){
2121 if(chargeOthers>0){
2122 return true;
2123 }
2124 else if(chargeOthers<0){
2125 return false;
2126 }
7099c89a 2127 }
fa8e4941 2128
2129 else{
28cc9dd5 2130 AliInfo("Error: Charge not lower nor higher as zero");
fa8e4941 2131 return false;
7099c89a 2132 }
fa8e4941 2133
28cc9dd5 2134 AliInfo("Error: Check values of Charge");
7099c89a 2135 return false;
7099c89a 2136}
fa8e4941 2137
2138