]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/UserTasks/AliAnalysisTaskPPJetSpectra.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
- Added special multiplicity estimator for pp
[u/mrichter/AliRoot.git] / PWGJE / UserTasks / AliAnalysisTaskPPJetSpectra.cxx
CommitLineData
3cfb058d 1#include "AliAnalysisTaskSE.h"
2#include "AliAnalysisHelperJetTasks.h"
3#include "THnSparse.h"
4#include "AliAnalysisManager.h"
5#include "AliInputEventHandler.h"
6#include "AliAODHandler.h"
7#include "AliAODInputHandler.h"
8#include "AliAODEvent.h"
9#include "AliESDEvent.h"
10#include "TAxis.h"
11#include "AliAODMCParticle.h"
12
13
14#include "AliAnalysisTaskPPJetSpectra.h"
15#include "AliAnalysisHelperJetTasks.h"
16
17ClassImp(AliAnalysisTaskPPJetSpectra)
18
19//---------------------------------------------------------------------------------------------------
20AliAnalysisTaskPPJetSpectra::AliAnalysisTaskPPJetSpectra()
21 :AliAnalysisTaskSE()
22 ,fOutputList(0)
23 ,fESD(0)
24 ,fAOD(0)
dc477024 25 ,fAODIn(0)
26 ,fAODOut(0)
27 ,fAODExt(0)
3cfb058d 28 ,fNonStdFile("")
29 ,fDebug(0)
30 ,fUseMC(kFALSE)
31 ,fEvtSelectionMask(0)
32 ,fEventClass(-1)
33 ,nVtxContCut(2)
34 ,fVtxZcut(10)
35 ,fVtxRcut(1)
36 ,nTrackFilter(272)
37 ,trackPtMin(0.15)
38 ,trackPtMax(100.)
39 ,trackEtaAbsMax(0.9)
40 ,jetPtMin(0.)
41 ,jetEtaCut(0.)
42 ,jetZmax(0.)
43 ,fhnEvent(0)
44 ,fhnTracks(0)
45 ,fhnMC(0)
46 ,fhnMC2(0)
47 ,fhnRecJetsNoCut(0)
48 ,fhnGenJetsNoCut(0)
49 ,fhnRecJetsCut(0)
50 ,fhnGenJetsCut(0)
51 ,fhnRecBckg(0)
52 ,fhnGenBckg(0)
53 ,fhnRecJetsTrackUEcor(0)
54 ,fhnGenJetsTrackUEcor(0)
55 ,fhnRecJetsBckgUEcor(0)
56 ,fhnGenJetsBckgUEcor(0)
57 ,fhnTrackUE(0)
58 ,fhnParticleUE(0)
59 ,fhnBckgRecUE(0)
60 ,fhnBckgGenUE(0)
61 ,fRecJetBranch("")
62 ,fGenJetBranch("")
63 ,fRecBckgBranch("")
64 ,fGenBckgBranch("")
65 ,fRecJetR(0.)
66 ,fGenJetR(0.)
67 ,fRecBckgR(0.)
68 ,fGenBckgR(0.)
69 ,fTrackType(0)
70 ,fParticleType(0)
71 ,fhnMatching(0)
72 ,fhnTrackCorrMatching(0)
73 ,fhnBckgCorrMatching(0)
74 ,fhnTrackUEanal(0)
75 ,kDoUEanalysis(kFALSE)
76 ,fRejectPileUp(1)
77{
78 if(fDebug) printf("%s: %d Constructor\n",(char*)__FILE__, __LINE__);
79 //DefineOutput(1, TList::Class());
80}
81
82//---------------------------------------------------------------------------------------------------
83AliAnalysisTaskPPJetSpectra::AliAnalysisTaskPPJetSpectra(const char* name):AliAnalysisTaskSE(name)
84 ,fOutputList(0)
85 ,fESD(0)
86 ,fAOD(0)
dc477024 87 ,fAODIn(0)
88 ,fAODOut(0)
89 ,fAODExt(0)
3cfb058d 90 ,fNonStdFile("")
91 ,fDebug(0)
92 ,fUseMC(kFALSE)
93 ,fEvtSelectionMask(0)
94 ,fEventClass(-1)
95 ,nVtxContCut(2)
96 ,fVtxZcut(10)
97 ,fVtxRcut(1)
98 ,nTrackFilter(272)
99 ,trackPtMin(0.15)
100 ,trackPtMax(100.)
101 ,trackEtaAbsMax(0.9)
102 ,jetPtMin(0.)
103 ,jetEtaCut(0.)
104 ,jetZmax(0.)
105 ,fhnEvent(0)
106 ,fhnTracks(0)
107 ,fhnMC(0)
108 ,fhnMC2(0)
109 ,fhnRecJetsNoCut(0)
110 ,fhnGenJetsNoCut(0)
111 ,fhnRecJetsCut(0)
112 ,fhnGenJetsCut(0)
113 ,fhnRecBckg(0)
114 ,fhnGenBckg(0)
115 ,fhnRecJetsTrackUEcor(0)
116 ,fhnGenJetsTrackUEcor(0)
117 ,fhnRecJetsBckgUEcor(0)
118 ,fhnGenJetsBckgUEcor(0)
119 ,fhnTrackUE(0)
120 ,fhnParticleUE(0)
121 ,fhnBckgRecUE(0)
122 ,fhnBckgGenUE(0)
123 ,fRecJetBranch("")
124 ,fGenJetBranch("")
125 ,fRecBckgBranch("")
126 ,fGenBckgBranch("")
127 ,fRecJetR(0.)
128 ,fGenJetR(0.)
129 ,fRecBckgR(0.)
130 ,fGenBckgR(0.)
131 ,fTrackType(0)
132 ,fParticleType(0)
133 ,fhnMatching(0)
134 ,fhnTrackCorrMatching(0)
135 ,fhnBckgCorrMatching(0)
136 ,fhnTrackUEanal(0)
137 ,kDoUEanalysis(kFALSE)
138 ,fRejectPileUp(1)
139{
140 if(fDebug) printf("%s: %d Constructor\n",(char*)__FILE__, __LINE__);
141 DefineOutput(1, TList::Class());
142}
143
144//---------------------------------------------------------------------------------------------------
145void AliAnalysisTaskPPJetSpectra::Init()
146{
147 if(fDebug) printf("%s: %d Initialize\n", (char*)__FILE__, __LINE__);
148 return;
149}
150
151//---------------------------------------------------------------------------------------------------
152Bool_t AliAnalysisTaskPPJetSpectra::Notify()
153{
154
155 if(fDebug) printf("%s: %d Notify\n", (char*)__FILE__, __LINE__);
156 return kTRUE;
157}
158
159//---------------------------------------------------------------------------------------------------
160void AliAnalysisTaskPPJetSpectra::UserCreateOutputObjects()
161{
162 fOutputList = new TList;
163 PostData(1,fOutputList);
164 fOutputList->SetOwner(1);
165
166 const Double_t binsPt[] = {0.,2.,4.,6.,8.,10.,12., 14.,16., 18.,20.,24.,28.,32.,38.,44.,50.,58.,66.,76.,86.,100.,114.,128.,150.,250.};
167 const Int_t nBinsPt = 25;
168
169 const Double_t trackPtBins[42] = {0,0.2,0.4,0.6,0.8,1.,1.2,1.4,1.6,1.8,2.,2.5,3.,3.5,4.,4.5,5.,5.5,6.,8.,10.,12.,14.,16.,18.,20.,25.,30.,35.,40.,45.,50.,55.,60.,65.,70.,75.,80.,85.,90.,95.,100};
170 const Int_t nTrackPtBins = 41;
171
172 if(fDebug) printf("%s: %d Creating output objects\n",(char*)__FILE__,__LINE__);
173 // triggered - centrality - z-vtxs - r-vtx - ncontrib - ntracks - IsPileUp?- njetsrec - njetsgen - accepted
174 const Int_t evtsDim = 7;
175 Int_t evtsBins[evtsDim] = { 2, 22, 40, 20, 20, 2, 2};
176 Double_t evtsMins[evtsDim] = {-0.5, -5., -20., 0., 0., -0.5, -0.5};
177 Double_t evtsMaxs[evtsDim] = { 1.5, 105., 20., 2., 40., 1.5, 1.5};
178 // pt - eta - phi - is272 - charge
179 const Int_t trksDim = 7;
180 Int_t trksBins[trksDim] = { nTrackPtBins, 10, 10, 2, 2, 2, 3};
181 Double_t trksMins[trksDim] = { 0., -1., 0., -0.5,-0.5,-0.5, -1.5};
182 Double_t trksMaxs[trksDim] = { 0., 1., 2*TMath::Pi(), 1.5, 1.5, 1.5, 1.5};
183 // pt - eta - phi - AreaCh - AreaN - ptLeadingTrack/ptJet - nTracks
184 const Int_t jetsDim = 7;
185 Int_t jetsBins[jetsDim] = { nBinsPt, 20, 10,3, 5, 20, 30};
186 Double_t jetsMins[jetsDim] = { 0., -1., 0., 0., 0., 0., 0.};
187 Double_t jetsMaxs[jetsDim] = { 0., 1., 2*TMath::Pi(), 3., 5., 1.,30.};
188 // pt - eta - phi - AreaCH - AreaN - UEdensity - ptLeadingTrack/ptJet - ptOrig
189 const Int_t jetsUEcorDim = 8;
190 Int_t jetsUEcorBins[jetsUEcorDim] = { nBinsPt, 10, 10, 3, 5, 100, 20, nBinsPt};
191 Double_t jetsUEcorMins[jetsUEcorDim] = { 0., -1., 0., 0., 0., 0., 0., 0.};
192 Double_t jetsUEcorMaxs[jetsUEcorDim] = { 0., 1., 2*TMath::Pi(), 3., 5., 10., 1., 0.};
193 // leading jet pt - eta - sumAllPt1 - sumAllPt2 - Rpar
194 const Int_t ueDim = 5;
195 Int_t ueBins[ueDim] = { nBinsPt, 20, 100, 100, 100};
196 Double_t ueMins[ueDim] = { 0., -1, 0., 0., 0.};
197 Double_t ueMaxs[ueDim] = { 0., 1, 10., 10., 10.};
198 // ptGen - ptRec - ptFraction - maxDist - bothWays - (1 - rec/gen) - dR - dEta - dPhi - isSamePtBin
199 const Int_t matchingDim = 10;
200 Int_t matchingBins[matchingDim] = {nBinsPt, nBinsPt, 20, 5, 2, 100, 10, 10, 10, 2 };
201 Double_t matchingMins[matchingDim] = { 0., 0., 0., 0.05, -0.5, 0, 0, 0, 0,-0.5};
202 Double_t matchingMaxs[matchingDim] = { 0., 0., 1., 0.55, 1.5, 1,0.5, 0.5,0.5, 1.5};
203
204 fhnEvent = new THnSparseF("fhnEvent", "fhnEvent;triggered;centrality;zVtx;rVtx;Ncontrib;IsPileUp;isAccepted", evtsDim, evtsBins, evtsMins, evtsMaxs);
205 fhnTracks = new THnSparseF("fhnTracks", "fhnTraks;pT;eta;phi;isFilterMask272;isMask16;isMask256;charge;", trksDim, trksBins, trksMins, trksMaxs);
206 fhnMC = new THnSparseF("fhnMC", "fhnTraks;pT;eta;phi;isFilterMask272;isMask16;isMask256;charge;", trksDim, trksBins, trksMins, trksMaxs);
207 fhnMC2 = new THnSparseF("fhnMC2", "fhnTraks;pT;eta;phi;isFilterMask272;isMask16;isMask256;charge;", trksDim, trksBins, trksMins, trksMaxs);
208
209 fhnRecJetsNoCut = new THnSparseF("fhnRecJetsNoCut", "fhnRecJetsNoCut;pT;eta;phi;AreaCH;AreaN;z;nTracks;", jetsDim, jetsBins, jetsMins, jetsMaxs);
210 fhnGenJetsNoCut = new THnSparseF("fhnGenJetsNoCut", "fhnGenJetsNoCut;pT;eta;phi;AreaCH;AreaN;z;nTracks;", jetsDim, jetsBins, jetsMins, jetsMaxs);
211 fhnRecJetsCut = new THnSparseF("fhnRecJetsCut", "fhnRecJetsCut;pT;eta;phi;AreaCH;AreaN;z;nTracks;", jetsDim, jetsBins, jetsMins, jetsMaxs);
212 fhnGenJetsCut = new THnSparseF("fhnGenJetsCut", "fhnGenJetsCut;pT;eta;phi;AreaCH;AreaN;z;nTracks;", jetsDim, jetsBins, jetsMins, jetsMaxs);
213 //UE sparses
214 fhnTrackUE = new THnSparseF("fhnTrackUE", "fhnTrackUE;pT-leading jet;eta;sum1;sum2;avgSum;", ueDim, ueBins, ueMins, ueMaxs);
215 fhnParticleUE = new THnSparseF("fhnParticleUE","fhnParticleUE;pT-leading jet;eta;sum1;sum2;avgSum;", ueDim, ueBins, ueMins, ueMaxs);
216 // pt - eta - phi - AreaCH - AreaN - UEdensity - ptLeadingTrack/ptJet - ptOrig
217 fhnRecJetsTrackUEcor = new THnSparseF("fhnRecJetsTrackUEcor", "fhnRecJetsTrackUEcor;pT;eta;phi;AreaCH;AreaN;UE;z;pTorig;", jetsUEcorDim, jetsUEcorBins, jetsUEcorMins, jetsUEcorMaxs);
218 fhnGenJetsTrackUEcor = new THnSparseF("fhnGenJetsTrackUEcor", "fhnGenJetsTrackUEcor;pT;eta;phi;AreaCH;AreaN;UE;z;pTorig;", jetsUEcorDim, jetsUEcorBins, jetsUEcorMins, jetsUEcorMaxs);
219 // ptGen - ptRec - ptFraction - maxDist - bothWays
220 fhnMatching = new THnSparseF("fhnMatching","fhnMatching;pTgen;pTrec;ptFraction;maxDist;bothWays;1-p_{T,rec}/p_{T,gen}; dR; dEta; dPhi;;isSamePtBin;", matchingDim, matchingBins, matchingMins, matchingMaxs);
221 fhnTrackCorrMatching = new THnSparseF("fhnTrackCorrMatching","fhnTrackCorrMatching;pTgen;pTrec;ptFraction;maxDist;bothWays;1-p_{T,rec}/p_{T,gen}; dR; dEta; dPhi;isSamePtBin;", matchingDim, matchingBins, matchingMins, matchingMaxs);
222
223 fhnTracks->SetBinEdges(0, trackPtBins);
224 fhnTrackUE->SetBinEdges(0, trackPtBins);
225 fhnParticleUE->SetBinEdges(0, trackPtBins);
226
227 fhnRecJetsNoCut->SetBinEdges(0, binsPt);
228 fhnGenJetsNoCut->SetBinEdges(0, binsPt);
229 fhnRecJetsCut->SetBinEdges(0, binsPt);
230 fhnGenJetsCut->SetBinEdges(0, binsPt);
231
232 fhnRecJetsTrackUEcor->SetBinEdges(0,binsPt);
233 fhnRecJetsTrackUEcor->SetBinEdges(7,binsPt);
234 fhnGenJetsTrackUEcor->SetBinEdges(0,binsPt);
235 fhnGenJetsTrackUEcor->SetBinEdges(7,binsPt);
236
237 fhnMatching->SetBinEdges(0,binsPt);
238 fhnMatching->SetBinEdges(1,binsPt);
239 fhnTrackCorrMatching->SetBinEdges(0,binsPt);
240 fhnTrackCorrMatching->SetBinEdges(1,binsPt);
241
242 fhnEvent->Sumw2();
243 fhnTracks->Sumw2();
244 fhnMC->Sumw2();
245 fhnMC2->Sumw2();
246 fhnRecJetsNoCut->Sumw2();
247 fhnGenJetsNoCut->Sumw2();
248 fhnRecJetsCut->Sumw2();
249 fhnGenJetsCut->Sumw2();
250 fhnRecJetsTrackUEcor->Sumw2();
251 fhnGenJetsTrackUEcor->Sumw2();
252 fhnTrackUE->Sumw2();
253 fhnParticleUE->Sumw2();
254 fhnMatching->Sumw2();
255 fhnTrackCorrMatching->Sumw2();
256
257 fOutputList->Add(fhnEvent);
258 fOutputList->Add(fhnTracks);
259 fOutputList->Add(fhnMC);
260 fOutputList->Add(fhnMC2);
261 fOutputList->Add(fhnRecJetsNoCut);
262 fOutputList->Add(fhnGenJetsNoCut);
263 fOutputList->Add(fhnRecJetsCut);
264 fOutputList->Add(fhnGenJetsCut);
265 fOutputList->Add(fhnRecJetsTrackUEcor);
266 fOutputList->Add(fhnGenJetsTrackUEcor);
267 fOutputList->Add(fhnTrackUE);
268 fOutputList->Add(fhnParticleUE);
269 fOutputList->Add(fhnMatching);
270 fOutputList->Add(fhnTrackCorrMatching);
271
272 const Int_t trackUEdim = 11;
273 //max track pt - ch_part_dens - ch_pt_dens - avg_pt_track - zero particles -tMAX_part_dens-tMIN_part_dens-tDIF_part_dens-tMAX_pt_dens-tMIN_pt_dens-tDIF_pt_dens
274 Int_t trackUEbins[trackUEdim] = { 100, 50, 50, 50, 2, 50, 50, 50, 50, 50, 50};
275 Double_t trackUEmins[trackUEdim] = { 0., 0., 0., 0.,-0.5, 0., 0., 0., 0., 0., 0.};
276 Double_t trackUEmaxs[trackUEdim] = {100., 5., 5., 5., 1.5, 5., 5., 5., 5., 5., 5.};
277
278 fhnTrackUEanal = new THnSparseF("fhnTrackUEanal","fhnTrackUEanal;PTmax;#rho^{ch.part}_{T};#rho^{pT}_{T};avg-pt/track;zero_particles;#rho^{ch.part}_{transMAX};#rho^{ch.part}_{transMIN};#rho^{ch.part}_{transDIF};#rho^{pT}_{transMAX};#rho^{pT}_{transMIN};#rho^{pT}_{transDIF}", trackUEdim, trackUEbins, trackUEmins, trackUEmaxs);
279 fhnTrackUEanal->Sumw2();
280 fOutputList->Add(fhnTrackUEanal);
281
282}
283
284//---------------------------------------------------------------------------------------------------
dc477024 285void AliAnalysisTaskPPJetSpectra::UserExec(Option_t */*option*/)
3cfb058d 286{
287 Double_t evtContainer[7];
288 Bool_t isEventSelected = EventSelection(evtContainer);
289
290 TList trackList;
291 TList particleList;
292 Int_t nTracks = GetListOfTracks(fTrackType, &trackList);
293 Int_t nParticles = GetListOfTracks(fParticleType, &particleList);
294
295 if(fDebug > 10 )
296 std::cout<< nTracks<< " " << nParticles<<std::endl;
297
298 TClonesArray *tcaRecJets = 0;
299 TClonesArray *tcaGenJets = 0;
300 TClonesArray *tcaRecBckg = 0;
301 TClonesArray *tcaGenBckg = 0;
302
303 Int_t rJ, gJ, rB, gB;
dc477024 304 if(fAODOut && !tcaRecJets && fRecJetBranch.Length() > 0) tcaRecJets = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fRecJetBranch.Data()));
305 if(fAODOut && !tcaGenJets && fGenJetBranch.Length() > 0) tcaGenJets = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fGenJetBranch.Data()));
306 if(fAODOut && !tcaRecBckg && fRecBckgBranch.Length() > 0) tcaRecBckg = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fRecBckgBranch.Data()));
307 if(fAODOut && !tcaGenBckg && fGenBckgBranch.Length() > 0) tcaGenBckg = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fGenBckgBranch.Data()));
308
309 if(fAODExt) {
310 if(!tcaRecJets && fRecJetBranch.Length() > 0) tcaRecJets = dynamic_cast<TClonesArray*>(fAODExt->GetAOD()->FindListObject(fRecJetBranch.Data()));
311 if(!tcaGenJets && fGenJetBranch.Length() > 0) tcaGenJets = dynamic_cast<TClonesArray*>(fAODExt->GetAOD()->FindListObject(fGenJetBranch.Data()));
312 if(!tcaRecBckg && fRecBckgBranch.Length() > 0) tcaRecBckg = dynamic_cast<TClonesArray*>(fAODExt->GetAOD()->FindListObject(fRecBckgBranch.Data()));
313 if(!tcaGenBckg && fGenBckgBranch.Length() > 0) tcaGenBckg = dynamic_cast<TClonesArray*>(fAODExt->GetAOD()->FindListObject(fGenBckgBranch.Data()));
314 }
315
316 if(fAODIn) {
317 if(!tcaRecJets && fRecJetBranch.Length() > 0) tcaRecJets = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fRecJetBranch.Data()));
318 if(!tcaGenJets && fGenJetBranch.Length() > 0) tcaGenJets = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fGenJetBranch.Data()));
319 if(!tcaRecBckg && fRecBckgBranch.Length() > 0) tcaRecBckg = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fRecBckgBranch.Data()));
320 if(!tcaGenBckg && fGenBckgBranch.Length() > 0) tcaGenBckg = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fGenBckgBranch.Data()));
321 }
3cfb058d 322
dc477024 323 if(!tcaRecJets) rJ = -1;
324 else rJ = tcaRecJets->GetEntries();
325 if(!tcaGenJets) gJ = -1;
326 else gJ = tcaGenJets->GetEntries();
327 if(!tcaRecBckg) rB = -1;
328 else rB = tcaRecBckg->GetEntries();
329 if(!tcaGenBckg) gB = -1;
330 else gB = tcaGenBckg->GetEntries();
3cfb058d 331
332 TList jetRecListNoCut;
333 TList jetGenListNoCut;
334 Int_t nRecJetsNoCut = (rJ>0?GetListOfJets(tcaRecJets, &jetRecListNoCut, kFALSE):0);
335 Int_t nGenJetsNoCut = (gJ>0?GetListOfJets(tcaGenJets, &jetGenListNoCut, kFALSE):0);
336
337 TList jetRecListCut;
338 TList jetGenListCut;
339 Int_t nRecJetsCut = (rJ>0?GetListOfJets(tcaRecJets, &jetRecListCut, kTRUE):0);
340 Int_t nGenJetsCut = (gJ>0?GetListOfJets(tcaGenJets, &jetGenListCut, kTRUE):0);
341
342 TList bckgRecListCut;
343 TList bckgGenListCut;
344 Int_t nRecBckgCut = (rB>0?GetListOfJets(tcaRecBckg, &bckgRecListCut, kTRUE):0);
345 Int_t nGenBckgCut = (gB>0?GetListOfJets(tcaGenBckg, &bckgGenListCut, kTRUE):0);
346
347 if(fDebug > 10) std::cout<<nRecJetsNoCut<<" "<<nGenJetsNoCut<<" "<<nRecJetsCut<< " "<< nGenJetsCut<<" "<<nRecBckgCut<<" "<<nGenBckgCut<<std::endl;
348
349 evtContainer[6] = (isEventSelected?1:0);
350 fhnEvent->Fill(evtContainer);
351
352 if(!isEventSelected) {
353 PostData(1,fOutputList);
354 return;
355 }
356
357 if(kDoUEanalysis) DoUEAnalysis(&trackList, (Double_t)0.15, (Double_t)0.9);
358
359 Double_t recUEdensity = GetUE(&jetRecListCut, &trackList, fRecJetR, fhnTrackUE);
360 Double_t genUEdensity = GetUE(&jetGenListCut, &particleList, fGenJetR, fhnParticleUE);
361
362 TList trackUEcorrJetRecListCut;
363 TList trackUEcorrJetGenListCut;
364
365 CorrectForUE(&jetRecListCut, recUEdensity, &trackUEcorrJetRecListCut, fhnRecJetsTrackUEcor);
366 CorrectForUE(&jetGenListCut, genUEdensity, &trackUEcorrJetGenListCut, fhnGenJetsTrackUEcor);
367
368 FillJetContainer(&jetRecListNoCut, fhnRecJetsNoCut);
369 FillJetContainer(&jetGenListNoCut, fhnGenJetsNoCut);
370 FillJetContainer(&jetRecListCut, fhnRecJetsCut);
371 FillJetContainer(&jetGenListCut, fhnGenJetsCut);
372
373 if(trackUEcorrJetRecListCut.GetEntries() != 0 && trackUEcorrJetGenListCut.GetEntries() != 0 ) {
374 MatchJets(kTRUE, &trackUEcorrJetRecListCut, &trackUEcorrJetGenListCut, 0.1, fhnTrackCorrMatching);
375 MatchJets(kTRUE, &trackUEcorrJetRecListCut, &trackUEcorrJetGenListCut, 0.3, fhnTrackCorrMatching);
376 MatchJets(kTRUE, &trackUEcorrJetRecListCut, &trackUEcorrJetGenListCut, 0.5, fhnTrackCorrMatching);
377 }
378 if(jetRecListCut.GetEntries() != 0 && jetGenListCut.GetEntries() != 0 ) {
379 MatchJets(kTRUE, &jetRecListCut, &jetGenListCut, 0.1, fhnMatching);
380 MatchJets(kTRUE, &jetRecListCut, &jetGenListCut, 0.3, fhnMatching);
381 MatchJets(kTRUE, &jetRecListCut, &jetGenListCut, 0.5, fhnMatching);
382 }
383
384
385 PostData(1,fOutputList);
386 return;
387}
388
389//---------------------------------------------------------------------------------------------------
390void AliAnalysisTaskPPJetSpectra::Terminate(Option_t *)
391{
392 if(fDebug) printf("%s: %d Terminating\n",(char*)__FILE__, __LINE__);
393 return;
394}
395
396//---------------------------------------------------------------------------------------------------
397void AliAnalysisTaskPPJetSpectra::SetRecJetBranch(TString s) {
398 fRecJetBranch = s;
399 if(s.Contains("KT02") || s.Contains("SISCONE02") || s.Contains("UA102")) fRecJetR = 0.2;
400 else if(s.Contains("KT03") || s.Contains("SISCONE03") || s.Contains("UA103")) fRecJetR = 0.3;
401 else if(s.Contains("KT04") || s.Contains("SISCONE04") || s.Contains("UA104")) fRecJetR = 0.4;
402 else if(s.Contains("KT05") || s.Contains("SISCONE05") || s.Contains("UA105")) fRecJetR = 0.5;
403 else if(s.Contains("KT06") || s.Contains("SISCONE06") || s.Contains("UA106")) fRecJetR = 0.6;
404 else fRecJetR = 0;
405 return;
406}
407
408//---------------------------------------------------------------------------------------------------
409void AliAnalysisTaskPPJetSpectra::SetGenJetBranch(TString s) {
410 fGenJetBranch = s;
411 if(s.Contains("KT02") || s.Contains("SISCONE02") || s.Contains("UA102")) fGenJetR = 0.2;
412 else if(s.Contains("KT03") || s.Contains("SISCONE03") || s.Contains("UA103")) fGenJetR = 0.3;
413 else if(s.Contains("KT04") || s.Contains("SISCONE04") || s.Contains("UA104")) fGenJetR = 0.4;
414 else if(s.Contains("KT05") || s.Contains("SISCONE05") || s.Contains("UA105")) fGenJetR = 0.5;
415 else if(s.Contains("KT06") || s.Contains("SISCONE06") || s.Contains("UA106")) fGenJetR = 0.6;
416 else fGenJetR = 0;
417 return;
418}
419
420//---------------------------------------------------------------------------------------------------
421void AliAnalysisTaskPPJetSpectra::SetRecBckgBranch(TString s) {
422 fRecBckgBranch = s;
423 if(s.Contains("KT02") || s.Contains("SISCONE02") || s.Contains("UA102")) fRecBckgR = 0.2;
424 else if(s.Contains("KT03") || s.Contains("SISCONE03") || s.Contains("UA103")) fRecBckgR = 0.3;
425 else if(s.Contains("KT04") || s.Contains("SISCONE04") || s.Contains("UA104")) fRecBckgR = 0.4;
426 else if(s.Contains("KT05") || s.Contains("SISCONE05") || s.Contains("UA105")) fRecBckgR = 0.5;
427 else if(s.Contains("KT06") || s.Contains("SISCONE06") || s.Contains("UA106")) fRecBckgR = 0.6;
428 else fRecBckgR = 0;
429 return;
430}
431
432//---------------------------------------------------------------------------------------------------
433void AliAnalysisTaskPPJetSpectra::SetGenBckgBranch(TString s) {
434 fGenBckgBranch = s;
435 if(s.Contains("KT02") || s.Contains("SISCONE02") || s.Contains("UA102")) fGenBckgR = 0.2;
436 else if(s.Contains("KT03") || s.Contains("SISCONE03") || s.Contains("UA103")) fGenBckgR = 0.3;
437 else if(s.Contains("KT04") || s.Contains("SISCONE04") || s.Contains("UA104")) fGenBckgR = 0.4;
438 else if(s.Contains("KT05") || s.Contains("SISCONE05") || s.Contains("UA105")) fGenBckgR = 0.5;
439 else if(s.Contains("KT06") || s.Contains("SISCONE06") || s.Contains("UA106")) fGenBckgR = 0.6;
440 else fGenBckgR = 0;
441 return;
442}
443
444//---------------------------------------------------------------------------------------------------
445void AliAnalysisTaskPPJetSpectra::SetVertexCuts(Int_t nCont, Double_t Zcut, Double_t Rcut) {
446 if(fDebug)
447 printf("%s: %d Setting Vertex Cuts\n",(char*)__FILE__,__LINE__);
448 nVtxContCut = nCont;
449 fVtxZcut = Zcut;
450 fVtxRcut = Rcut;
451 return;
452}
453
454//---------------------------------------------------------------------------------------------------
455void AliAnalysisTaskPPJetSpectra::SetTrackCuts(Double_t ptMin, Double_t ptMax, Double_t etaMax) {
456 if(fDebug > 1)
457 printf("%s: %d Track cuts set\n",(char*)__FILE__,__LINE__);
458 trackPtMin = ptMin;
459 trackPtMax = ptMax;
460 trackEtaAbsMax = etaMax;
461 return;
462}
463
464//---------------------------------------------------------------------------------------------------
465void AliAnalysisTaskPPJetSpectra::SetJetCuts(Double_t ptMin, Double_t eta, Double_t z) {
466 jetPtMin = ptMin;
467 jetEtaCut = eta;
468 jetZmax = z;
469 return;
470}
471
472//---------------------------------------------------------------------------------------------------
473Bool_t AliAnalysisTaskPPJetSpectra::EventSelection(Double_t evtContainer[6]) {
474 if(fDebug > 1)
475 printf("%s: %d UserExec analysing event %d.\n", (char*)__FILE__, __LINE__, (Int_t)fEntry);
476
477 // Trigger selection
478 AliInputEventHandler* inputHandler = (AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
479 if(fDebug > 2 )
480 printf("IsEventSelected: %d eventSelectionMask: %d\n", (Int_t)inputHandler->IsEventSelected(), (Int_t)fEvtSelectionMask);
481
482 Bool_t isTriggerSelected = inputHandler->IsEventSelected() & fEvtSelectionMask;
483 evtContainer[0] = (Int_t)isTriggerSelected;
484
dc477024 485 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
3cfb058d 486 if(!fESD && fDebug > 2)
487 printf("%s: %d No ESD event found\n",(char*)__FILE__, __LINE__);
488
dc477024 489 if(!fESD) fAODIn = dynamic_cast<AliAODEvent*>(InputEvent());
490 fAODOut = dynamic_cast<AliAODEvent*>(AODEvent());
3cfb058d 491
dc477024 492 if(!fESD)fAOD = fAODIn;
493 else fAOD = fAODOut;
3cfb058d 494
dc477024 495 if(fNonStdFile.Length()!=0) {
496 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
497 fAODExt = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
3cfb058d 498 }
dc477024 499
500 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
3cfb058d 501 Float_t centrality = -1;
502 if(fEventClass > 0)
503 {
dc477024 504 if(handler->InheritsFrom("AliAODHandler")) centrality = fAODIn->GetHeader()->GetCentrality();
3cfb058d 505 else if(fESD) centrality = fESD->GetCentrality()->GetCentralityPercentile("V0M");
506 else centrality = AliAnalysisHelperJetTasks::EventClass();
507 }
508 evtContainer[1] = centrality;
509
510 // Event selection
511 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
512 if(!primVtx) {
513 evtContainer[2] = 1e10;
514 evtContainer[3] = 1e10;
515 evtContainer[4] = 1e10;
516 return kFALSE;
517 } else {
518 primVtx->Print();
519 evtContainer[2] = primVtx->GetZ();
520 evtContainer[3] = TMath::Sqrt(primVtx->GetX()*primVtx->GetX() + primVtx->GetY()*primVtx->GetY() );
521 evtContainer[4] = primVtx->GetNContributors();
522 }
523
524 if(fRejectPileUp && AliAnalysisHelperJetTasks::IsPileUp()){
525 if (fDebug > 1) Printf("%s:%d SPD pileup: event REJECTED...", (char*)__FILE__,__LINE__);
526 evtContainer[5] = 0.;
527 return kFALSE;
528 }
529 else
530 evtContainer[5] = 1.;
531
532 if(evtContainer[0] == 0)
533 {
534 if(fDebug > 2)
535 printf("%s: %d Event rejected: Trigger\n",(char*)__FILE__, __LINE__);
536 return kFALSE;
537 }
538 if(evtContainer[1] >= 0 && evtContainer[1] < 10*((Int_t)fEventClass/10) && evtContainer[1] > 10*((Int_t)fEventClass/10 + 1) )
539 {
540 if(fDebug > 2)
541 printf("%s: %d Event rejected: Centrality\n",(char*)__FILE__, __LINE__);
542 return kFALSE;
543 }
544 if(TMath::Abs(evtContainer[2]) > fVtxZcut)
545 {
546 if(fDebug > 2)
547 printf("%s: %d Event rejected: Vertex position (z)\n",(char*)__FILE__, __LINE__);
548 return kFALSE;
549 }
550 if(evtContainer[3] > fVtxRcut)
551 {
552 if(fDebug > 2)
553 printf("%s: %d Event rejected: Vertex position (x-y)\n",(char*)__FILE__,__LINE__);
554 return kFALSE;
555 }
556 if(evtContainer[4] < nVtxContCut)
557 {
558 if(fDebug > 2)
559 printf("%s: %d Event rejected: Vertex contributors \n",(char*)__FILE__, __LINE__);
560 return kFALSE;
561 }
562
563 return kTRUE;
564}
565
566//---------------------------------------------------------------------------------------------------
567Int_t AliAnalysisTaskPPJetSpectra::GetListOfTracks(Int_t trackType, TList *trackList) {
568 if(fDebug)
569 printf("%s: %d Getting Tracks\n",(char*)__FILE__,__LINE__);
570 if(trackType == AliAnalysisTaskPPJetSpectra::kNone) {
571 if(fDebug) printf("%s: %d No track type selected", (char*)__FILE__,__LINE__);
572 return 0;
573 }
574
575 if(1 == trackType) {
576 for(Int_t i = 0; i < fAOD->GetNumberOfTracks(); i++) {
577 Bool_t isOK = kTRUE;
578 Bool_t isOK16 = kTRUE;
579 Bool_t isOK256 = kTRUE;
580 AliAODTrack* track = (AliAODTrack*)fAOD->GetTrack(i);
581 if(!track) continue;
582 if(nTrackFilter>0 && !track->TestFilterBit(nTrackFilter) ) isOK = kFALSE;
583 if(nTrackFilter>0 && !track->TestFilterBit(16) ) isOK16 = kFALSE;
584 if(nTrackFilter>0 && !track->TestFilterBit(256) ) isOK256 = kFALSE;
585 if(track->Pt() < trackPtMin || track->Pt() > trackPtMax) isOK = kFALSE;
586 if(trackEtaAbsMax < TMath::Abs(track->Eta())) isOK = kFALSE;
587 Double_t tmpContainer[] = {track->Pt(), track->Eta(), track->Phi(), (Double_t)isOK, (Double_t)isOK16, (Double_t)isOK256, (Double_t)track->Charge()/3};
588 if(isOK) trackList->Add(track);
589 fhnTracks->Fill(tmpContainer);
590 }
591 }
592 else {
593 TClonesArray *tca=dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
dc477024 594 if(!tca) {
595 if(fDebug > 2) printf("no branch %s\n", AliAODMCParticle::StdBranchName());
596 }
3cfb058d 597 for(Int_t i = 0; i < tca->GetEntries(); i++) {
598 AliAODMCParticle *particle = (AliAODMCParticle*)tca->At(i);
599 if(!particle) continue;
600 if(!particle->IsPhysicalPrimary()) continue;
601 Double_t charge = particle->Charge();
602 if( charge > 0 ) charge = 1;
603 if( charge < 0 ) charge = -1;
604 if( (trackType == 3) && (charge == 0) ) continue;
605 trackList->Add(particle);
606 Double_t tmpContainer[] = {particle->Pt(), particle->Eta(), particle->Phi(), -999, -999, -999, charge };
607 if(trackType == 2) fhnMC->Fill(tmpContainer);
608 else fhnMC2->Fill(tmpContainer);
609 }
610 }
611 trackList->Sort();
612 return trackList->GetEntries();
613}
614
3cfb058d 615//---------------------------------------------------------------------------------------------------
616Int_t AliAnalysisTaskPPJetSpectra::GetListOfJets(TClonesArray* jca, TList* jetList, Bool_t useCuts) {
617 if(fDebug)
618 printf("%s: %d Getting list of jets\n",(char*)__FILE__,__LINE__);
619
620 Double_t ptCut = (useCuts && jetPtMin>0 ?jetPtMin :0.1);
621 Double_t etaCut = (useCuts && jetEtaCut>0?jetEtaCut:0.9);
622 Double_t zCut = (useCuts && jetZmax>0 ?jetZmax :1);
623
624 if(!jca){
625 return -2;
626 }
627
628 for(Int_t i =0 ; i < jca->GetEntries(); i++) {
629 AliAODJet* jet = (AliAODJet*)jca->At(i);
630 if(!jet) continue;
631 Double_t zLeading = -1;
632// if(jet->GetPtLeading() <= 0) {
633// TRefArray* jetRef = (TRefArray*)jet->GetRefTracks();
634// for(int j = 0; j < jetRef->GetEntries(); j++) {
635// AliVTrack *track = (AliVTrack*)jet->GetTrack(j);
636// if(!track) continue;
637// cout<<track->Pt()<<" "<<zLeading<<endl;
638// if(track->Pt() > zLeading) zLeading = track->Pt();
639// cout<<__LINE__<<endl;
640// }
641// if(jet->Pt() > 0) zLeading /= jet->Pt();
642// else zLeading = -1;
643// }
644 //cuty
645 if(ptCut > 0 && ptCut > jet->Pt() ) continue;
646 if(etaCut > 0 && etaCut < TMath::Abs(jet->Eta()) ) continue;
647 if(zCut > 0 && zCut < zLeading && zLeading > 0 ) continue;
648
649 jetList->Add(jet);
650 }
651
652 jetList->Sort();
653 return jetList->GetEntries();
654}
655
656//---------------------------------------------------------------------------------------------------
657void AliAnalysisTaskPPJetSpectra::FillJetContainer(TList* jets, THnSparseF* container) {
658 if(!jets || !container) return;
659
660 for(Int_t i =0 ; i < jets->GetEntries(); i++) {
661 AliAODJet* jet = (AliAODJet*)jets->At(i);
662 if(!jet) continue;
663 // pt - eta - phi - AreaCH - AreaN - ptLeadingTrack/ptJet - ptOrig
664 Double_t tmpContainer[] = {jet->Pt(), jet->Eta(), jet->Phi(), jet->EffectiveAreaCharged(), jet->EffectiveAreaNeutral(), jet->GetPtLeading()/jet->Pt(), (Double_t)((TRefArray*)jet->GetRefTracks())->GetEntries() };
665 container->Fill(tmpContainer);
666 }
667 return;
668}
669
670//---------------------------------------------------------------------------------------------------
671Double_t AliAnalysisTaskPPJetSpectra::GetUE(TList* jets, TList* tracks, Double_t Rpar,THnSparseF* thnCont) {
672 if(!jets || !tracks) {
673 if(fDebug) printf("%s: %d No jets or tracks\n",(char*)__FILE__,__LINE__);
674 return 0;
675 }
676
677 //find leading jet
678 AliAODJet *jet = 0;
679 for(Int_t iTmp = 0; iTmp < jets->GetEntries(); iTmp++) {
680 AliAODJet* tmp = (AliAODJet*)jets->At(iTmp);
681 if(!tmp) continue;
682 if(!jet) jet = tmp;
683 else
684 if( jet->Pt() < tmp->Pt() ) jet = tmp;
685 }
686 if(!jet) return 0;
687
688 //get track in perpendicular cone
689 //TVector3 for jet and perpendicular cones with same eta axis
690 TVector3 j(jet->Px(), jet->Py(), jet->Pz());
691 TVector3 p1(j);
692 TVector3 p2(j);
693
694 p1.RotateZ(TMath::Pi()/2.);
695 p2.RotateZ(-TMath::Pi()/2.);
696
697 Double_t sumAllPt1 = 0;
698 Double_t sumAllPt2 = 0;
699
700 for(int i = 0; i < tracks->GetEntries(); i++) {
701 AliVParticle* tr = (AliVParticle*)tracks->At(i);
702 if(!tr) continue;
703 TVector3 v(tr->Px(), tr->Py(), tr->Pz());
704 Double_t dR1 = v.DrEtaPhi(p1);
705 Double_t dR2 = v.DrEtaPhi(p2);
706
707 if(dR1 < Rpar) sumAllPt1+=v.Pt();
708 if(dR2 < Rpar) sumAllPt2+=v.Pt();
709 }
710
711 if(Rpar != 0) {
712 sumAllPt1/=(TMath::Pi()*Rpar*Rpar);
713 sumAllPt2/=(TMath::Pi()*Rpar*Rpar);
714 }
715
716 Double_t container[] = {j.Pt(), j.Eta(), sumAllPt1, sumAllPt2, (sumAllPt1+sumAllPt2)/2.};
717 thnCont->Fill(container);
718
719 return (sumAllPt1+sumAllPt2)/2.;
720
721}
722
723//---------------------------------------------------------------------------------------------------
724Double_t AliAnalysisTaskPPJetSpectra::GetBckgUE(TList* jets, Double_t Rpar, Bool_t isSkipped, THnSparseF* cont) {
725 if(!jets) return 0;
726
727 Int_t leading01 = -1;
728 Int_t leading02 = -1;
729 Double_t leading01pt = 0;
730 Double_t leading02pt = 0;
731
732 if(!isSkipped) {
733 for(Int_t i = 0; i < jets->GetEntries(); i++) {
734 AliAODJet* tmp = (AliAODJet*)jets->At(i);
735 if(!tmp) continue;
736 if(tmp->Pt() > leading01pt) {
737 leading02 = leading01;
738 leading01 = i;
739 leading02pt = leading01pt;
740 leading01pt = tmp->Pt();
741 } else if(tmp->Pt() > leading02) {
742 leading02 = i;
743 leading02pt = tmp->Pt();
744 }
745 }
746 }
747
748 Int_t subLeading01 = -1;
749 Int_t subLeading02 = -1;
750 Double_t subLeading01pt = 0;
751 Double_t subLeading02pt = 0;
752 Double_t subLeading01area = 0;
753 Double_t subLeading02area = 0;
754
755 Double_t avgAll = 0;
756 Double_t areaAll = 0;
757 //search for subleading jets assuming leading 2 (skip02) jets are skipped
758 if(jets->GetEntries() > 0) {
759 for(Int_t i = 0; i < jets->GetEntries(); i++) {
760 if(i == leading01 || i == leading02) continue;
761 AliAODJet* tmp = (AliAODJet*)jets->At(i);
762 if(!tmp) continue;
763 if( subLeading01pt < tmp->Pt() ) {
764 subLeading02 = subLeading01;
765 subLeading01 = i;
766 subLeading02pt = subLeading01pt;
767 subLeading01pt = tmp->Pt();
768 subLeading02area = subLeading01area;
769 subLeading01area = tmp->EffectiveAreaCharged();
770 } else if( subLeading02pt < tmp->Pt() ) {
771 subLeading02 = i;
772 subLeading02pt = tmp->Pt();
773 subLeading02area = tmp->EffectiveAreaCharged();
774 }
775 }
776
777 //do avg jet pt
778 Int_t nJets = 0;
779 for(int i = 0; i < jets->GetEntries();i++) {
780 if(i == leading01 || i == leading02) continue;
781 AliAODJet* tmp = (AliAODJet*)jets->At(i);
782 if(!tmp) continue;
783 nJets++;
784 avgAll+= tmp->Pt();
785 areaAll += tmp->EffectiveAreaCharged();
786 }
787 }
788
789
790 if(avgAll!=0) avgAll/=areaAll;
791 if(subLeading01pt!=0) subLeading01pt/=subLeading01area;
792 if(subLeading02pt!=0)subLeading02pt/=subLeading02area;
793 Double_t tmpContainer[]={subLeading01pt,subLeading02pt,avgAll,Rpar*10};
794
795 cont->Fill(tmpContainer);
796
797 return avgAll;
798}
799
800//---------------------------------------------------------------------------------------------------
801Int_t AliAnalysisTaskPPJetSpectra::CorrectForUE(TList* jets, Double_t UE, TList* newList, THnSparseF* container) {
802 if(!jets) return 0;
803
804 for(Int_t i = 0; i < jets->GetEntries(); i++) {
805 AliAODJet* tmpJet = (AliAODJet*)jets->At(i);
806 if(!tmpJet) continue;
807
808 AliAODJet* jet = (AliAODJet*)tmpJet->Clone(Form("%s%s",tmpJet->GetName(),"_clone"));
809 Double_t corrUE = UE*(jet->EffectiveAreaCharged() + jet->EffectiveAreaNeutral() );
810
811 if(jet->Pt() > corrUE) jet->SetPtEtaPhiM(jet->Pt()- corrUE, jet->Eta(), jet->Phi(), jet->M());
812 else jet->SetPtEtaPhiM(0, jet->Eta(), jet->Phi(), jet->M());
813 newList->Add(jet);
814
815 // pt - eta - phi - AreaCH - AreaN - UEdensity - ptLeadingTrack/ptJet - ptOrig
816 Double_t cont[8]={jet->Pt(), jet->Eta(), jet->Phi(), jet->EffectiveAreaCharged(),jet->EffectiveAreaNeutral(), corrUE, (jet->Pt() > 0?jet->GetPtLeading() / jet->Pt():-1),tmpJet->Pt()};
817 container->Fill(cont);
818 }
819 newList->Sort();
820 return newList->GetEntries();
821}
822
823//---------------------------------------------------------------------------------------------------
824void AliAnalysisTaskPPJetSpectra::MatchJets(Bool_t doBothWay, TList* recList, TList* genList, Float_t maxDist, THnSparseF* container) {
825 if(!genList || !recList) return;
826 Int_t mode = 1;
827
828 const Int_t kGenJets = genList->GetEntries();
829 const Int_t kRecJets = recList->GetEntries();
830
831 TArrayI iMatchIndex(kGenJets);
832 TArrayF fPtFraction(kGenJets);
833
834 TArrayI aGenIndex(kRecJets);
835 TArrayI aRecIndex(kGenJets);
836
837 Double_t genPt = -1;
838 Double_t recPt = -1;
839 Double_t ptFraction = -1;
840
841 //-- Closest jets to generated and checked for E-fraction
842 if(!doBothWay) {
843 AliAnalysisHelperJetTasks::GetJetMatching(genList, kGenJets, recList, kRecJets, iMatchIndex, fPtFraction, fDebug, maxDist, mode);
844 for(Int_t ig = 0; ig < kGenJets; ig++) {
845 //get gen-jet
846 AliAODJet* genJet = (AliAODJet*)genList->At(ig);
847 if(!genJet) continue;
848 genPt = genJet->Pt();
849
850 //get rec-jet
851 Int_t ir = iMatchIndex[ig];
852 if(ir<0||ir>=recList->GetEntries()) continue;
853 AliAODJet* recJet = (AliAODJet*)recList->At(ir);
854 if(!recJet) continue;
855 recPt = recJet->Pt();
856
857 ptFraction = fPtFraction[ig];
858
859 Double_t dR = recJet->DeltaR(genJet);
860 Double_t dEta = recJet->Eta() - genJet->Eta();
861 Double_t dPhi2 = dR*dR - dEta*dEta;
862
863 Bool_t isSame = kFALSE;
864 if(CheckPtBin(genPt) == CheckPtBin(recPt)) isSame = kTRUE;
865 // ptGen - ptRec - ptFraction - maxDist - bothWays - (1 - rec/gen) - dR - dEta - dPhi
866 Double_t cont[] = {genPt,recPt,ptFraction,maxDist,(Float_t)doBothWay, (genPt != 0 ? (1. - recPt/genPt):1), dR, dEta, (dPhi2>=0?TMath::Sqrt(dPhi2):-1), (Double_t)isSame};
867 container->Fill(cont);
868 }
869 }
870 //-- Closest rec jet to gen and vice-versa & check for E-fraction
871 else {
872 AliAnalysisHelperJetTasks::GetClosestJets(genList, kGenJets, recList, kRecJets, aGenIndex, aRecIndex,fDebug, maxDist);
873 for(Int_t ig = 0; ig < kGenJets; ig++) {
874 //get gen-jet
875 AliAODJet* genJet = (AliAODJet*)genList->At(ig);
876 if(!genJet) continue;
877 genPt = genJet->Pt();
878
879 //get rec-jet
880 Int_t ir = aRecIndex[ig];
881 if(ir>=0&&ir<recList->GetEntries()) {
882 AliAODJet *recJet = (AliAODJet*)recList->At(ir);
883 if(!recJet) continue;
884 recPt = recJet->Pt();
885 Double_t dR = recJet->DeltaR(genJet);
886 Double_t dEta = recJet->Eta() - genJet->Eta();
887 Double_t dPhi2 = dR*dR - dEta*dEta;
888
889 ptFraction = AliAnalysisHelperJetTasks::GetFractionOfJet(recJet,genJet,mode);
890
891 Bool_t isSame = kFALSE;
892 if(CheckPtBin(genPt) == CheckPtBin(recPt)) isSame = kTRUE;
893 // ptGen - ptRec - ptFraction - maxDist - bothWays - (1 - rec/gen) - dR - dEta - dPhi
894 Double_t cont[] = {genPt,recPt,ptFraction,maxDist,(Float_t)doBothWay, (genPt != 0 ? (1. - recPt/genPt):1), dR, dEta, (dPhi2>=0?TMath::Sqrt(dPhi2):-1), (Double_t)isSame };
895 container->Fill(cont);
896 }
897 }
898 }
899
900 return;
901}
902
903//---------------------------------------------------------------------------------------------------
904Int_t AliAnalysisTaskPPJetSpectra::CheckPtBin(Double_t pT) {
905 Double_t bins[] = {0.,2.,4.,6.,8.,10.,12., 14.,16., 18.,20.,24.,28.,32.,38.,44.,50.,58.,66.,76.,86.,100.,114.,128.,150.,250.};
906 Double_t nBins = 25;
907
908 for(Int_t i = 0; i <= nBins; i++) {
909 if(pT < bins[i] ) return i;
910 }
911 return nBins+1;
912}
913
914//---------------------------------------------------------------------------------------------------
915void AliAnalysisTaskPPJetSpectra::DoUEAnalysis(TList* tracks, Double_t PTcut, Double_t ETAcut) {
916 Double_t PTmax = 0;
917 int iMax = -1;
918 for(int i = 0; i < tracks->GetEntries(); i++) {
919 AliVParticle *part = (AliVParticle*)tracks->At(i);
920 if(!part) continue;
921
922 if(part->Pt() < PTcut) continue;
923 if(TMath::Abs(part->Eta()) > ETAcut ) continue;
924
925 if(PTmax < part->Pt()) {
926 PTmax = part->Pt();
927 iMax = i;
928 }
929 }
930
931 if(PTmax < trackPtMin) return;
932
933 Double_t ch_part_density = 0;
934 Double_t ch_part_pt_density = 0;
935 Double_t avg_track_pt = 0;
936 Double_t zero_ch_part = 0;
937 Double_t transMAX_part_density = 0;
938 Double_t transMIN_part_density = 0;
939 Double_t transMAX_pt_density = 0;
940 Double_t transMIN_pt_density = 0;
941 Double_t transDIF_pt_density = 0;
942 Double_t transDIF_part_density = 0;
943
944 Double_t A1_part_density = 0;
945 Double_t A2_part_density = 0;
946 Double_t A1_pt_density = 0;
947 Double_t A2_pt_density = 0;
948
949 AliVParticle *maxParticle = (AliVParticle*)tracks->At(iMax);
950
951 Double_t PHImax = maxParticle->Phi();
952
953 for(int i = 0; i < tracks->GetEntries(); i++) {
954 AliVParticle* part = (AliVParticle*)tracks->At(i);
955 if(!part) continue;
956 if(part->Pt() < PTcut) continue;
957
958 Bool_t isArea1 = kFALSE;
959 Bool_t isArea2 = kFALSE;
960
961 Double_t dPhi = part->Phi() - PHImax;
962 if(dPhi > TMath::Pi()) dPhi -= TMath::Pi();
963 if(dPhi < -TMath::Pi()) dPhi +=TMath::Pi();
964
965 if(TMath::Abs(dPhi) > TMath::Pi()/3 && TMath::Abs(dPhi) < 2*TMath::Pi()/3 ) {
966 if(dPhi>0) isArea1 = kTRUE;
967 else isArea2 = kTRUE;
968 }
969
970 if(!isArea1 && !isArea2) continue;
971
972 ch_part_density += 1;
973 ch_part_pt_density += part->Pt();
974 avg_track_pt += part->Pt();
975
976 if(isArea1) {
977 A1_part_density += 1;
978 A1_pt_density += part->Pt();
979 } else if(isArea2) {
980 A2_part_density += 1;
981 A2_pt_density += part->Pt();
982 }
983 }
984
985
986 if(ch_part_density > 0) avg_track_pt/=ch_part_density;
987 else zero_ch_part = 1;
988 ch_part_density*=3./(4.*TMath::Pi()*ETAcut);
989 ch_part_pt_density*= 3./(4.*TMath::Pi()*ETAcut);
990
991 A1_part_density*=3./(4.*TMath::Pi()*ETAcut);
992 A1_pt_density*= 3./(4.*TMath::Pi()*ETAcut);
993 A2_part_density*=3./(4.*TMath::Pi()*ETAcut);
994 A2_pt_density*= 3./(4.*TMath::Pi()*ETAcut);
995
996 if(A1_part_density < A2_part_density) {
997 transMAX_part_density = A2_part_density;
998 transMIN_part_density = A1_part_density;
999 } else {
1000 transMAX_part_density = A1_part_density;
1001 transMIN_part_density = A2_part_density;
1002 }
1003
1004 if(A1_pt_density < A2_pt_density) {
1005 transMAX_pt_density = A2_pt_density;
1006 transMIN_pt_density = A1_pt_density;
1007 } else {
1008 transMAX_pt_density = A1_pt_density;
1009 transMIN_pt_density = A2_pt_density;
1010 }
1011
1012 transDIF_pt_density = transMAX_pt_density - transMIN_pt_density;
1013 transDIF_part_density = transMAX_part_density - transMIN_part_density;
1014
1015 Double_t container[] = {PTmax, ch_part_density, ch_part_pt_density, avg_track_pt, zero_ch_part,
1016 transMAX_part_density, transMIN_part_density, transDIF_part_density,
1017 transMAX_pt_density, transMIN_pt_density, transDIF_pt_density};
1018
1019 fhnTrackUEanal->Fill(container);
1020}