]> 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
Properly apply range at axis limits (Diego)
[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)
25 ,fMC(0)
26 ,fAODJets(0)
27 ,fAODExtension(0)
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)
87 ,fMC(0)
88 ,fAODJets(0)
89 ,fAODExtension(0)
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//---------------------------------------------------------------------------------------------------
285void AliAnalysisTaskPPJetSpectra::UserExec(Option_t *option)
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;
304
305 if(fRecJetBranch.Length() > 0) {
306 tcaRecJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fRecJetBranch.Data()));
307 if(!tcaRecJets) tcaRecJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fRecJetBranch.Data()));
308 if(!tcaRecJets) tcaRecJets = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fRecJetBranch.Data()));
309 if(!tcaRecJets) tcaRecJets = dynamic_cast<TClonesArray*>(fAOD->GetList()->FindObject(fRecJetBranch.Data()));
310 if(fAODExtension&&!tcaRecJets) tcaRecJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fRecJetBranch.Data()));
311 if(tcaRecJets) rJ = tcaRecJets->GetEntries();
312 else rJ = -1;
313 } else rJ = -1;
314
315 if(fGenJetBranch.Length() > 0) {
316 tcaGenJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fGenJetBranch.Data()));
317 if(!tcaGenJets) tcaGenJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fGenJetBranch.Data()));
318 if(!tcaGenJets) tcaGenJets = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fGenJetBranch.Data()));
319 if(!tcaGenJets) tcaGenJets = dynamic_cast<TClonesArray*>(fAOD->GetList()->FindObject(fGenJetBranch.Data()));
320 if(fAODExtension&&!tcaGenJets) tcaGenJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fGenJetBranch.Data()));
321 if(tcaGenJets) gJ = tcaGenJets->GetEntries();
322 else gJ = -1;
323 } else gJ = -1;
324
325 if(fRecBckgBranch.Length() > 0) {
326 tcaRecBckg = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fRecBckgBranch.Data()));
327 if(!tcaRecBckg) tcaRecBckg = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fRecBckgBranch.Data()));
328 if(!tcaRecBckg) tcaRecBckg = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fRecBckgBranch.Data()));
329 if(!tcaRecBckg) tcaRecBckg = dynamic_cast<TClonesArray*>(fAOD->GetList()->FindObject(fRecBckgBranch.Data()));
330 if(fAODExtension&&!tcaRecBckg) tcaRecBckg = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fRecBckgBranch.Data()));
331 rB = tcaRecBckg->GetEntries();
332 } else rB = -1;
333
334 if(fGenBckgBranch.Length() > 0) {
335 tcaGenBckg = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fGenBckgBranch.Data()));
336 if(!tcaGenBckg) tcaGenBckg = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fGenBckgBranch.Data()));
337 if(!tcaGenBckg) tcaGenBckg = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fGenBckgBranch.Data()));
338 if(!tcaGenBckg) tcaGenBckg = dynamic_cast<TClonesArray*>(fAOD->GetList()->FindObject(fGenBckgBranch.Data()));
339 if(fAODExtension&&!tcaGenBckg) tcaGenBckg = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fGenBckgBranch.Data()));
340 gB = tcaGenBckg->GetEntries();
341 } else gB = -1;
342
343 TList jetRecListNoCut;
344 TList jetGenListNoCut;
345 Int_t nRecJetsNoCut = (rJ>0?GetListOfJets(tcaRecJets, &jetRecListNoCut, kFALSE):0);
346 Int_t nGenJetsNoCut = (gJ>0?GetListOfJets(tcaGenJets, &jetGenListNoCut, kFALSE):0);
347
348 TList jetRecListCut;
349 TList jetGenListCut;
350 Int_t nRecJetsCut = (rJ>0?GetListOfJets(tcaRecJets, &jetRecListCut, kTRUE):0);
351 Int_t nGenJetsCut = (gJ>0?GetListOfJets(tcaGenJets, &jetGenListCut, kTRUE):0);
352
353 TList bckgRecListCut;
354 TList bckgGenListCut;
355 Int_t nRecBckgCut = (rB>0?GetListOfJets(tcaRecBckg, &bckgRecListCut, kTRUE):0);
356 Int_t nGenBckgCut = (gB>0?GetListOfJets(tcaGenBckg, &bckgGenListCut, kTRUE):0);
357
358 if(fDebug > 10) std::cout<<nRecJetsNoCut<<" "<<nGenJetsNoCut<<" "<<nRecJetsCut<< " "<< nGenJetsCut<<" "<<nRecBckgCut<<" "<<nGenBckgCut<<std::endl;
359
360 evtContainer[6] = (isEventSelected?1:0);
361 fhnEvent->Fill(evtContainer);
362
363 if(!isEventSelected) {
364 PostData(1,fOutputList);
365 return;
366 }
367
368 if(kDoUEanalysis) DoUEAnalysis(&trackList, (Double_t)0.15, (Double_t)0.9);
369
370 Double_t recUEdensity = GetUE(&jetRecListCut, &trackList, fRecJetR, fhnTrackUE);
371 Double_t genUEdensity = GetUE(&jetGenListCut, &particleList, fGenJetR, fhnParticleUE);
372
373 TList trackUEcorrJetRecListCut;
374 TList trackUEcorrJetGenListCut;
375
376 CorrectForUE(&jetRecListCut, recUEdensity, &trackUEcorrJetRecListCut, fhnRecJetsTrackUEcor);
377 CorrectForUE(&jetGenListCut, genUEdensity, &trackUEcorrJetGenListCut, fhnGenJetsTrackUEcor);
378
379 FillJetContainer(&jetRecListNoCut, fhnRecJetsNoCut);
380 FillJetContainer(&jetGenListNoCut, fhnGenJetsNoCut);
381 FillJetContainer(&jetRecListCut, fhnRecJetsCut);
382 FillJetContainer(&jetGenListCut, fhnGenJetsCut);
383
384 if(trackUEcorrJetRecListCut.GetEntries() != 0 && trackUEcorrJetGenListCut.GetEntries() != 0 ) {
385 MatchJets(kTRUE, &trackUEcorrJetRecListCut, &trackUEcorrJetGenListCut, 0.1, fhnTrackCorrMatching);
386 MatchJets(kTRUE, &trackUEcorrJetRecListCut, &trackUEcorrJetGenListCut, 0.3, fhnTrackCorrMatching);
387 MatchJets(kTRUE, &trackUEcorrJetRecListCut, &trackUEcorrJetGenListCut, 0.5, fhnTrackCorrMatching);
388 }
389 if(jetRecListCut.GetEntries() != 0 && jetGenListCut.GetEntries() != 0 ) {
390 MatchJets(kTRUE, &jetRecListCut, &jetGenListCut, 0.1, fhnMatching);
391 MatchJets(kTRUE, &jetRecListCut, &jetGenListCut, 0.3, fhnMatching);
392 MatchJets(kTRUE, &jetRecListCut, &jetGenListCut, 0.5, fhnMatching);
393 }
394
395
396 PostData(1,fOutputList);
397 return;
398}
399
400//---------------------------------------------------------------------------------------------------
401void AliAnalysisTaskPPJetSpectra::Terminate(Option_t *)
402{
403 if(fDebug) printf("%s: %d Terminating\n",(char*)__FILE__, __LINE__);
404 return;
405}
406
407//---------------------------------------------------------------------------------------------------
408void AliAnalysisTaskPPJetSpectra::SetRecJetBranch(TString s) {
409 fRecJetBranch = s;
410 if(s.Contains("KT02") || s.Contains("SISCONE02") || s.Contains("UA102")) fRecJetR = 0.2;
411 else if(s.Contains("KT03") || s.Contains("SISCONE03") || s.Contains("UA103")) fRecJetR = 0.3;
412 else if(s.Contains("KT04") || s.Contains("SISCONE04") || s.Contains("UA104")) fRecJetR = 0.4;
413 else if(s.Contains("KT05") || s.Contains("SISCONE05") || s.Contains("UA105")) fRecJetR = 0.5;
414 else if(s.Contains("KT06") || s.Contains("SISCONE06") || s.Contains("UA106")) fRecJetR = 0.6;
415 else fRecJetR = 0;
416 return;
417}
418
419//---------------------------------------------------------------------------------------------------
420void AliAnalysisTaskPPJetSpectra::SetGenJetBranch(TString s) {
421 fGenJetBranch = s;
422 if(s.Contains("KT02") || s.Contains("SISCONE02") || s.Contains("UA102")) fGenJetR = 0.2;
423 else if(s.Contains("KT03") || s.Contains("SISCONE03") || s.Contains("UA103")) fGenJetR = 0.3;
424 else if(s.Contains("KT04") || s.Contains("SISCONE04") || s.Contains("UA104")) fGenJetR = 0.4;
425 else if(s.Contains("KT05") || s.Contains("SISCONE05") || s.Contains("UA105")) fGenJetR = 0.5;
426 else if(s.Contains("KT06") || s.Contains("SISCONE06") || s.Contains("UA106")) fGenJetR = 0.6;
427 else fGenJetR = 0;
428 return;
429}
430
431//---------------------------------------------------------------------------------------------------
432void AliAnalysisTaskPPJetSpectra::SetRecBckgBranch(TString s) {
433 fRecBckgBranch = s;
434 if(s.Contains("KT02") || s.Contains("SISCONE02") || s.Contains("UA102")) fRecBckgR = 0.2;
435 else if(s.Contains("KT03") || s.Contains("SISCONE03") || s.Contains("UA103")) fRecBckgR = 0.3;
436 else if(s.Contains("KT04") || s.Contains("SISCONE04") || s.Contains("UA104")) fRecBckgR = 0.4;
437 else if(s.Contains("KT05") || s.Contains("SISCONE05") || s.Contains("UA105")) fRecBckgR = 0.5;
438 else if(s.Contains("KT06") || s.Contains("SISCONE06") || s.Contains("UA106")) fRecBckgR = 0.6;
439 else fRecBckgR = 0;
440 return;
441}
442
443//---------------------------------------------------------------------------------------------------
444void AliAnalysisTaskPPJetSpectra::SetGenBckgBranch(TString s) {
445 fGenBckgBranch = s;
446 if(s.Contains("KT02") || s.Contains("SISCONE02") || s.Contains("UA102")) fGenBckgR = 0.2;
447 else if(s.Contains("KT03") || s.Contains("SISCONE03") || s.Contains("UA103")) fGenBckgR = 0.3;
448 else if(s.Contains("KT04") || s.Contains("SISCONE04") || s.Contains("UA104")) fGenBckgR = 0.4;
449 else if(s.Contains("KT05") || s.Contains("SISCONE05") || s.Contains("UA105")) fGenBckgR = 0.5;
450 else if(s.Contains("KT06") || s.Contains("SISCONE06") || s.Contains("UA106")) fGenBckgR = 0.6;
451 else fGenBckgR = 0;
452 return;
453}
454
455//---------------------------------------------------------------------------------------------------
456void AliAnalysisTaskPPJetSpectra::SetVertexCuts(Int_t nCont, Double_t Zcut, Double_t Rcut) {
457 if(fDebug)
458 printf("%s: %d Setting Vertex Cuts\n",(char*)__FILE__,__LINE__);
459 nVtxContCut = nCont;
460 fVtxZcut = Zcut;
461 fVtxRcut = Rcut;
462 return;
463}
464
465//---------------------------------------------------------------------------------------------------
466void AliAnalysisTaskPPJetSpectra::SetTrackCuts(Double_t ptMin, Double_t ptMax, Double_t etaMax) {
467 if(fDebug > 1)
468 printf("%s: %d Track cuts set\n",(char*)__FILE__,__LINE__);
469 trackPtMin = ptMin;
470 trackPtMax = ptMax;
471 trackEtaAbsMax = etaMax;
472 return;
473}
474
475//---------------------------------------------------------------------------------------------------
476void AliAnalysisTaskPPJetSpectra::SetJetCuts(Double_t ptMin, Double_t eta, Double_t z) {
477 jetPtMin = ptMin;
478 jetEtaCut = eta;
479 jetZmax = z;
480 return;
481}
482
483//---------------------------------------------------------------------------------------------------
484Bool_t AliAnalysisTaskPPJetSpectra::EventSelection(Double_t evtContainer[6]) {
485 if(fDebug > 1)
486 printf("%s: %d UserExec analysing event %d.\n", (char*)__FILE__, __LINE__, (Int_t)fEntry);
487
488 // Trigger selection
489 AliInputEventHandler* inputHandler = (AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
490 if(fDebug > 2 )
491 printf("IsEventSelected: %d eventSelectionMask: %d\n", (Int_t)inputHandler->IsEventSelected(), (Int_t)fEvtSelectionMask);
492
493 Bool_t isTriggerSelected = inputHandler->IsEventSelected() & fEvtSelectionMask;
494 evtContainer[0] = (Int_t)isTriggerSelected;
495
496 fESD = (AliESDEvent*)InputEvent();
497 if(!fESD && fDebug > 2)
498 printf("%s: %d No ESD event found\n",(char*)__FILE__, __LINE__);
499
500 fMC = MCEvent();
501 if(!fMC && fDebug > 2)
502 printf("%s: %d No MC event found\n", (char*)__FILE__, __LINE__);
503
504
505 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
506 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
507 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
508 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
509 fAODJets = ((AliAODHandler*)handler)->GetAOD();
510 if (fDebug > 1)
511 Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
512 }
513 else {
514 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
515 if( handler && handler->InheritsFrom("AliAODHandler") ) {
516 fAOD = ((AliAODHandler*)handler)->GetAOD();
517 fAODJets = fAOD;
518 if (fDebug > 1)
519 Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
520 }
521 }
522
523 if(!fAODJets)
524 {
525 if(fDebug > 1)
526 printf("%s: %d No AOD found\n",(char*)__FILE__,__LINE__);
527 return kFALSE;
528 }
529 handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
530 Float_t centrality = -1;
531 if(fEventClass > 0)
532 {
533 if(handler->InheritsFrom("AliAODHandler")) centrality = fAOD->GetHeader()->GetCentrality();
534 else if(fESD) centrality = fESD->GetCentrality()->GetCentralityPercentile("V0M");
535 else centrality = AliAnalysisHelperJetTasks::EventClass();
536 }
537 evtContainer[1] = centrality;
538
539 // Event selection
540 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
541 if(!primVtx) {
542 evtContainer[2] = 1e10;
543 evtContainer[3] = 1e10;
544 evtContainer[4] = 1e10;
545 return kFALSE;
546 } else {
547 primVtx->Print();
548 evtContainer[2] = primVtx->GetZ();
549 evtContainer[3] = TMath::Sqrt(primVtx->GetX()*primVtx->GetX() + primVtx->GetY()*primVtx->GetY() );
550 evtContainer[4] = primVtx->GetNContributors();
551 }
552
553 if(fRejectPileUp && AliAnalysisHelperJetTasks::IsPileUp()){
554 if (fDebug > 1) Printf("%s:%d SPD pileup: event REJECTED...", (char*)__FILE__,__LINE__);
555 evtContainer[5] = 0.;
556 return kFALSE;
557 }
558 else
559 evtContainer[5] = 1.;
560
561 if(evtContainer[0] == 0)
562 {
563 if(fDebug > 2)
564 printf("%s: %d Event rejected: Trigger\n",(char*)__FILE__, __LINE__);
565 return kFALSE;
566 }
567 if(evtContainer[1] >= 0 && evtContainer[1] < 10*((Int_t)fEventClass/10) && evtContainer[1] > 10*((Int_t)fEventClass/10 + 1) )
568 {
569 if(fDebug > 2)
570 printf("%s: %d Event rejected: Centrality\n",(char*)__FILE__, __LINE__);
571 return kFALSE;
572 }
573 if(TMath::Abs(evtContainer[2]) > fVtxZcut)
574 {
575 if(fDebug > 2)
576 printf("%s: %d Event rejected: Vertex position (z)\n",(char*)__FILE__, __LINE__);
577 return kFALSE;
578 }
579 if(evtContainer[3] > fVtxRcut)
580 {
581 if(fDebug > 2)
582 printf("%s: %d Event rejected: Vertex position (x-y)\n",(char*)__FILE__,__LINE__);
583 return kFALSE;
584 }
585 if(evtContainer[4] < nVtxContCut)
586 {
587 if(fDebug > 2)
588 printf("%s: %d Event rejected: Vertex contributors \n",(char*)__FILE__, __LINE__);
589 return kFALSE;
590 }
591
592 return kTRUE;
593}
594
595//---------------------------------------------------------------------------------------------------
596Int_t AliAnalysisTaskPPJetSpectra::GetListOfTracks(Int_t trackType, TList *trackList) {
597 if(fDebug)
598 printf("%s: %d Getting Tracks\n",(char*)__FILE__,__LINE__);
599 if(trackType == AliAnalysisTaskPPJetSpectra::kNone) {
600 if(fDebug) printf("%s: %d No track type selected", (char*)__FILE__,__LINE__);
601 return 0;
602 }
603
604 if(1 == trackType) {
605 for(Int_t i = 0; i < fAOD->GetNumberOfTracks(); i++) {
606 Bool_t isOK = kTRUE;
607 Bool_t isOK16 = kTRUE;
608 Bool_t isOK256 = kTRUE;
609 AliAODTrack* track = (AliAODTrack*)fAOD->GetTrack(i);
610 if(!track) continue;
611 if(nTrackFilter>0 && !track->TestFilterBit(nTrackFilter) ) isOK = kFALSE;
612 if(nTrackFilter>0 && !track->TestFilterBit(16) ) isOK16 = kFALSE;
613 if(nTrackFilter>0 && !track->TestFilterBit(256) ) isOK256 = kFALSE;
614 if(track->Pt() < trackPtMin || track->Pt() > trackPtMax) isOK = kFALSE;
615 if(trackEtaAbsMax < TMath::Abs(track->Eta())) isOK = kFALSE;
616 Double_t tmpContainer[] = {track->Pt(), track->Eta(), track->Phi(), (Double_t)isOK, (Double_t)isOK16, (Double_t)isOK256, (Double_t)track->Charge()/3};
617 if(isOK) trackList->Add(track);
618 fhnTracks->Fill(tmpContainer);
619 }
620 }
621 else {
622 TClonesArray *tca=dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
623 for(Int_t i = 0; i < tca->GetEntries(); i++) {
624 AliAODMCParticle *particle = (AliAODMCParticle*)tca->At(i);
625 if(!particle) continue;
626 if(!particle->IsPhysicalPrimary()) continue;
627 Double_t charge = particle->Charge();
628 if( charge > 0 ) charge = 1;
629 if( charge < 0 ) charge = -1;
630 if( (trackType == 3) && (charge == 0) ) continue;
631 trackList->Add(particle);
632 Double_t tmpContainer[] = {particle->Pt(), particle->Eta(), particle->Phi(), -999, -999, -999, charge };
633 if(trackType == 2) fhnMC->Fill(tmpContainer);
634 else fhnMC2->Fill(tmpContainer);
635 }
636 }
637 trackList->Sort();
638 return trackList->GetEntries();
639}
640
641//---------------------------------------------------------------------------------------------------
642Int_t AliAnalysisTaskPPJetSpectra::GetTCAJets(TClonesArray* jca, TString branch) {
643 if( branch.Length() == 0) {
644 if(fDebug) printf("%s: %d No branch \"%s\" selected\n", (char*)__FILE__,__LINE__, branch.Data() );
645 return 0;
646 }
647
648 jca = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(branch.Data()));
649 if(!jca) jca = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(branch.Data()));
650 if(!jca) jca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(branch.Data()));
651 if(!jca) jca = dynamic_cast<TClonesArray*>(fAOD->GetList()->FindObject(branch.Data()));
652 if(fAODExtension&&!jca) jca = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(branch.Data()));
653 return jca->GetEntries();
654
655}
656
657//---------------------------------------------------------------------------------------------------
658Int_t AliAnalysisTaskPPJetSpectra::GetListOfJets(TClonesArray* jca, TList* jetList, Bool_t useCuts) {
659 if(fDebug)
660 printf("%s: %d Getting list of jets\n",(char*)__FILE__,__LINE__);
661
662 Double_t ptCut = (useCuts && jetPtMin>0 ?jetPtMin :0.1);
663 Double_t etaCut = (useCuts && jetEtaCut>0?jetEtaCut:0.9);
664 Double_t zCut = (useCuts && jetZmax>0 ?jetZmax :1);
665
666 if(!jca){
667 return -2;
668 }
669
670 for(Int_t i =0 ; i < jca->GetEntries(); i++) {
671 AliAODJet* jet = (AliAODJet*)jca->At(i);
672 if(!jet) continue;
673 Double_t zLeading = -1;
674// if(jet->GetPtLeading() <= 0) {
675// TRefArray* jetRef = (TRefArray*)jet->GetRefTracks();
676// for(int j = 0; j < jetRef->GetEntries(); j++) {
677// AliVTrack *track = (AliVTrack*)jet->GetTrack(j);
678// if(!track) continue;
679// cout<<track->Pt()<<" "<<zLeading<<endl;
680// if(track->Pt() > zLeading) zLeading = track->Pt();
681// cout<<__LINE__<<endl;
682// }
683// if(jet->Pt() > 0) zLeading /= jet->Pt();
684// else zLeading = -1;
685// }
686 //cuty
687 if(ptCut > 0 && ptCut > jet->Pt() ) continue;
688 if(etaCut > 0 && etaCut < TMath::Abs(jet->Eta()) ) continue;
689 if(zCut > 0 && zCut < zLeading && zLeading > 0 ) continue;
690
691 jetList->Add(jet);
692 }
693
694 jetList->Sort();
695 return jetList->GetEntries();
696}
697
698//---------------------------------------------------------------------------------------------------
699void AliAnalysisTaskPPJetSpectra::FillJetContainer(TList* jets, THnSparseF* container) {
700 if(!jets || !container) return;
701
702 for(Int_t i =0 ; i < jets->GetEntries(); i++) {
703 AliAODJet* jet = (AliAODJet*)jets->At(i);
704 if(!jet) continue;
705 // pt - eta - phi - AreaCH - AreaN - ptLeadingTrack/ptJet - ptOrig
706 Double_t tmpContainer[] = {jet->Pt(), jet->Eta(), jet->Phi(), jet->EffectiveAreaCharged(), jet->EffectiveAreaNeutral(), jet->GetPtLeading()/jet->Pt(), (Double_t)((TRefArray*)jet->GetRefTracks())->GetEntries() };
707 container->Fill(tmpContainer);
708 }
709 return;
710}
711
712//---------------------------------------------------------------------------------------------------
713Double_t AliAnalysisTaskPPJetSpectra::GetUE(TList* jets, TList* tracks, Double_t Rpar,THnSparseF* thnCont) {
714 if(!jets || !tracks) {
715 if(fDebug) printf("%s: %d No jets or tracks\n",(char*)__FILE__,__LINE__);
716 return 0;
717 }
718
719 //find leading jet
720 AliAODJet *jet = 0;
721 for(Int_t iTmp = 0; iTmp < jets->GetEntries(); iTmp++) {
722 AliAODJet* tmp = (AliAODJet*)jets->At(iTmp);
723 if(!tmp) continue;
724 if(!jet) jet = tmp;
725 else
726 if( jet->Pt() < tmp->Pt() ) jet = tmp;
727 }
728 if(!jet) return 0;
729
730 //get track in perpendicular cone
731 //TVector3 for jet and perpendicular cones with same eta axis
732 TVector3 j(jet->Px(), jet->Py(), jet->Pz());
733 TVector3 p1(j);
734 TVector3 p2(j);
735
736 p1.RotateZ(TMath::Pi()/2.);
737 p2.RotateZ(-TMath::Pi()/2.);
738
739 Double_t sumAllPt1 = 0;
740 Double_t sumAllPt2 = 0;
741
742 for(int i = 0; i < tracks->GetEntries(); i++) {
743 AliVParticle* tr = (AliVParticle*)tracks->At(i);
744 if(!tr) continue;
745 TVector3 v(tr->Px(), tr->Py(), tr->Pz());
746 Double_t dR1 = v.DrEtaPhi(p1);
747 Double_t dR2 = v.DrEtaPhi(p2);
748
749 if(dR1 < Rpar) sumAllPt1+=v.Pt();
750 if(dR2 < Rpar) sumAllPt2+=v.Pt();
751 }
752
753 if(Rpar != 0) {
754 sumAllPt1/=(TMath::Pi()*Rpar*Rpar);
755 sumAllPt2/=(TMath::Pi()*Rpar*Rpar);
756 }
757
758 Double_t container[] = {j.Pt(), j.Eta(), sumAllPt1, sumAllPt2, (sumAllPt1+sumAllPt2)/2.};
759 thnCont->Fill(container);
760
761 return (sumAllPt1+sumAllPt2)/2.;
762
763}
764
765//---------------------------------------------------------------------------------------------------
766Double_t AliAnalysisTaskPPJetSpectra::GetBckgUE(TList* jets, Double_t Rpar, Bool_t isSkipped, THnSparseF* cont) {
767 if(!jets) return 0;
768
769 Int_t leading01 = -1;
770 Int_t leading02 = -1;
771 Double_t leading01pt = 0;
772 Double_t leading02pt = 0;
773
774 if(!isSkipped) {
775 for(Int_t i = 0; i < jets->GetEntries(); i++) {
776 AliAODJet* tmp = (AliAODJet*)jets->At(i);
777 if(!tmp) continue;
778 if(tmp->Pt() > leading01pt) {
779 leading02 = leading01;
780 leading01 = i;
781 leading02pt = leading01pt;
782 leading01pt = tmp->Pt();
783 } else if(tmp->Pt() > leading02) {
784 leading02 = i;
785 leading02pt = tmp->Pt();
786 }
787 }
788 }
789
790 Int_t subLeading01 = -1;
791 Int_t subLeading02 = -1;
792 Double_t subLeading01pt = 0;
793 Double_t subLeading02pt = 0;
794 Double_t subLeading01area = 0;
795 Double_t subLeading02area = 0;
796
797 Double_t avgAll = 0;
798 Double_t areaAll = 0;
799 //search for subleading jets assuming leading 2 (skip02) jets are skipped
800 if(jets->GetEntries() > 0) {
801 for(Int_t i = 0; i < jets->GetEntries(); i++) {
802 if(i == leading01 || i == leading02) continue;
803 AliAODJet* tmp = (AliAODJet*)jets->At(i);
804 if(!tmp) continue;
805 if( subLeading01pt < tmp->Pt() ) {
806 subLeading02 = subLeading01;
807 subLeading01 = i;
808 subLeading02pt = subLeading01pt;
809 subLeading01pt = tmp->Pt();
810 subLeading02area = subLeading01area;
811 subLeading01area = tmp->EffectiveAreaCharged();
812 } else if( subLeading02pt < tmp->Pt() ) {
813 subLeading02 = i;
814 subLeading02pt = tmp->Pt();
815 subLeading02area = tmp->EffectiveAreaCharged();
816 }
817 }
818
819 //do avg jet pt
820 Int_t nJets = 0;
821 for(int i = 0; i < jets->GetEntries();i++) {
822 if(i == leading01 || i == leading02) continue;
823 AliAODJet* tmp = (AliAODJet*)jets->At(i);
824 if(!tmp) continue;
825 nJets++;
826 avgAll+= tmp->Pt();
827 areaAll += tmp->EffectiveAreaCharged();
828 }
829 }
830
831
832 if(avgAll!=0) avgAll/=areaAll;
833 if(subLeading01pt!=0) subLeading01pt/=subLeading01area;
834 if(subLeading02pt!=0)subLeading02pt/=subLeading02area;
835 Double_t tmpContainer[]={subLeading01pt,subLeading02pt,avgAll,Rpar*10};
836
837 cont->Fill(tmpContainer);
838
839 return avgAll;
840}
841
842//---------------------------------------------------------------------------------------------------
843Int_t AliAnalysisTaskPPJetSpectra::CorrectForUE(TList* jets, Double_t UE, TList* newList, THnSparseF* container) {
844 if(!jets) return 0;
845
846 for(Int_t i = 0; i < jets->GetEntries(); i++) {
847 AliAODJet* tmpJet = (AliAODJet*)jets->At(i);
848 if(!tmpJet) continue;
849
850 AliAODJet* jet = (AliAODJet*)tmpJet->Clone(Form("%s%s",tmpJet->GetName(),"_clone"));
851 Double_t corrUE = UE*(jet->EffectiveAreaCharged() + jet->EffectiveAreaNeutral() );
852
853 if(jet->Pt() > corrUE) jet->SetPtEtaPhiM(jet->Pt()- corrUE, jet->Eta(), jet->Phi(), jet->M());
854 else jet->SetPtEtaPhiM(0, jet->Eta(), jet->Phi(), jet->M());
855 newList->Add(jet);
856
857 // pt - eta - phi - AreaCH - AreaN - UEdensity - ptLeadingTrack/ptJet - ptOrig
858 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()};
859 container->Fill(cont);
860 }
861 newList->Sort();
862 return newList->GetEntries();
863}
864
865//---------------------------------------------------------------------------------------------------
866void AliAnalysisTaskPPJetSpectra::MatchJets(Bool_t doBothWay, TList* recList, TList* genList, Float_t maxDist, THnSparseF* container) {
867 if(!genList || !recList) return;
868 Int_t mode = 1;
869
870 const Int_t kGenJets = genList->GetEntries();
871 const Int_t kRecJets = recList->GetEntries();
872
873 TArrayI iMatchIndex(kGenJets);
874 TArrayF fPtFraction(kGenJets);
875
876 TArrayI aGenIndex(kRecJets);
877 TArrayI aRecIndex(kGenJets);
878
879 Double_t genPt = -1;
880 Double_t recPt = -1;
881 Double_t ptFraction = -1;
882
883 //-- Closest jets to generated and checked for E-fraction
884 if(!doBothWay) {
885 AliAnalysisHelperJetTasks::GetJetMatching(genList, kGenJets, recList, kRecJets, iMatchIndex, fPtFraction, fDebug, maxDist, mode);
886 for(Int_t ig = 0; ig < kGenJets; ig++) {
887 //get gen-jet
888 AliAODJet* genJet = (AliAODJet*)genList->At(ig);
889 if(!genJet) continue;
890 genPt = genJet->Pt();
891
892 //get rec-jet
893 Int_t ir = iMatchIndex[ig];
894 if(ir<0||ir>=recList->GetEntries()) continue;
895 AliAODJet* recJet = (AliAODJet*)recList->At(ir);
896 if(!recJet) continue;
897 recPt = recJet->Pt();
898
899 ptFraction = fPtFraction[ig];
900
901 Double_t dR = recJet->DeltaR(genJet);
902 Double_t dEta = recJet->Eta() - genJet->Eta();
903 Double_t dPhi2 = dR*dR - dEta*dEta;
904
905 Bool_t isSame = kFALSE;
906 if(CheckPtBin(genPt) == CheckPtBin(recPt)) isSame = kTRUE;
907 // ptGen - ptRec - ptFraction - maxDist - bothWays - (1 - rec/gen) - dR - dEta - dPhi
908 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};
909 container->Fill(cont);
910 }
911 }
912 //-- Closest rec jet to gen and vice-versa & check for E-fraction
913 else {
914 AliAnalysisHelperJetTasks::GetClosestJets(genList, kGenJets, recList, kRecJets, aGenIndex, aRecIndex,fDebug, maxDist);
915 for(Int_t ig = 0; ig < kGenJets; ig++) {
916 //get gen-jet
917 AliAODJet* genJet = (AliAODJet*)genList->At(ig);
918 if(!genJet) continue;
919 genPt = genJet->Pt();
920
921 //get rec-jet
922 Int_t ir = aRecIndex[ig];
923 if(ir>=0&&ir<recList->GetEntries()) {
924 AliAODJet *recJet = (AliAODJet*)recList->At(ir);
925 if(!recJet) continue;
926 recPt = recJet->Pt();
927 Double_t dR = recJet->DeltaR(genJet);
928 Double_t dEta = recJet->Eta() - genJet->Eta();
929 Double_t dPhi2 = dR*dR - dEta*dEta;
930
931 ptFraction = AliAnalysisHelperJetTasks::GetFractionOfJet(recJet,genJet,mode);
932
933 Bool_t isSame = kFALSE;
934 if(CheckPtBin(genPt) == CheckPtBin(recPt)) isSame = kTRUE;
935 // ptGen - ptRec - ptFraction - maxDist - bothWays - (1 - rec/gen) - dR - dEta - dPhi
936 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 };
937 container->Fill(cont);
938 }
939 }
940 }
941
942 return;
943}
944
945//---------------------------------------------------------------------------------------------------
946Int_t AliAnalysisTaskPPJetSpectra::CheckPtBin(Double_t pT) {
947 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.};
948 Double_t nBins = 25;
949
950 for(Int_t i = 0; i <= nBins; i++) {
951 if(pT < bins[i] ) return i;
952 }
953 return nBins+1;
954}
955
956//---------------------------------------------------------------------------------------------------
957void AliAnalysisTaskPPJetSpectra::DoUEAnalysis(TList* tracks, Double_t PTcut, Double_t ETAcut) {
958 Double_t PTmax = 0;
959 int iMax = -1;
960 for(int i = 0; i < tracks->GetEntries(); i++) {
961 AliVParticle *part = (AliVParticle*)tracks->At(i);
962 if(!part) continue;
963
964 if(part->Pt() < PTcut) continue;
965 if(TMath::Abs(part->Eta()) > ETAcut ) continue;
966
967 if(PTmax < part->Pt()) {
968 PTmax = part->Pt();
969 iMax = i;
970 }
971 }
972
973 if(PTmax < trackPtMin) return;
974
975 Double_t ch_part_density = 0;
976 Double_t ch_part_pt_density = 0;
977 Double_t avg_track_pt = 0;
978 Double_t zero_ch_part = 0;
979 Double_t transMAX_part_density = 0;
980 Double_t transMIN_part_density = 0;
981 Double_t transMAX_pt_density = 0;
982 Double_t transMIN_pt_density = 0;
983 Double_t transDIF_pt_density = 0;
984 Double_t transDIF_part_density = 0;
985
986 Double_t A1_part_density = 0;
987 Double_t A2_part_density = 0;
988 Double_t A1_pt_density = 0;
989 Double_t A2_pt_density = 0;
990
991 AliVParticle *maxParticle = (AliVParticle*)tracks->At(iMax);
992
993 Double_t PHImax = maxParticle->Phi();
994
995 for(int i = 0; i < tracks->GetEntries(); i++) {
996 AliVParticle* part = (AliVParticle*)tracks->At(i);
997 if(!part) continue;
998 if(part->Pt() < PTcut) continue;
999
1000 Bool_t isArea1 = kFALSE;
1001 Bool_t isArea2 = kFALSE;
1002
1003 Double_t dPhi = part->Phi() - PHImax;
1004 if(dPhi > TMath::Pi()) dPhi -= TMath::Pi();
1005 if(dPhi < -TMath::Pi()) dPhi +=TMath::Pi();
1006
1007 if(TMath::Abs(dPhi) > TMath::Pi()/3 && TMath::Abs(dPhi) < 2*TMath::Pi()/3 ) {
1008 if(dPhi>0) isArea1 = kTRUE;
1009 else isArea2 = kTRUE;
1010 }
1011
1012 if(!isArea1 && !isArea2) continue;
1013
1014 ch_part_density += 1;
1015 ch_part_pt_density += part->Pt();
1016 avg_track_pt += part->Pt();
1017
1018 if(isArea1) {
1019 A1_part_density += 1;
1020 A1_pt_density += part->Pt();
1021 } else if(isArea2) {
1022 A2_part_density += 1;
1023 A2_pt_density += part->Pt();
1024 }
1025 }
1026
1027
1028 if(ch_part_density > 0) avg_track_pt/=ch_part_density;
1029 else zero_ch_part = 1;
1030 ch_part_density*=3./(4.*TMath::Pi()*ETAcut);
1031 ch_part_pt_density*= 3./(4.*TMath::Pi()*ETAcut);
1032
1033 A1_part_density*=3./(4.*TMath::Pi()*ETAcut);
1034 A1_pt_density*= 3./(4.*TMath::Pi()*ETAcut);
1035 A2_part_density*=3./(4.*TMath::Pi()*ETAcut);
1036 A2_pt_density*= 3./(4.*TMath::Pi()*ETAcut);
1037
1038 if(A1_part_density < A2_part_density) {
1039 transMAX_part_density = A2_part_density;
1040 transMIN_part_density = A1_part_density;
1041 } else {
1042 transMAX_part_density = A1_part_density;
1043 transMIN_part_density = A2_part_density;
1044 }
1045
1046 if(A1_pt_density < A2_pt_density) {
1047 transMAX_pt_density = A2_pt_density;
1048 transMIN_pt_density = A1_pt_density;
1049 } else {
1050 transMAX_pt_density = A1_pt_density;
1051 transMIN_pt_density = A2_pt_density;
1052 }
1053
1054 transDIF_pt_density = transMAX_pt_density - transMIN_pt_density;
1055 transDIF_part_density = transMAX_part_density - transMIN_part_density;
1056
1057 Double_t container[] = {PTmax, ch_part_density, ch_part_pt_density, avg_track_pt, zero_ch_part,
1058 transMAX_part_density, transMIN_part_density, transDIF_part_density,
1059 transMAX_pt_density, transMIN_pt_density, transDIF_pt_density};
1060
1061 fhnTrackUEanal->Fill(container);
1062}