7 #include <TClonesArray.h>
9 #include <TLorentzVector.h>
11 #include <TParameter.h>
13 #include "AliAnalysisTask.h"
14 #include "AliAnalysisManager.h"
15 #include "AliEmcalJet.h"
16 #include "AliVCluster.h"
18 #include "AliESDEvent.h"
19 #include "AliAODEvent.h"
20 #include "AliESDInputHandler.h"
21 #include "AliCentrality.h"
23 #include "AliAnalysisTaskEmcalJetSpectra.h"
25 ClassImp(AliAnalysisTaskEmcalJetSpectra)
27 //________________________________________________________________________
28 AliAnalysisTaskEmcalJetSpectra::AliAnalysisTaskEmcalJetSpectra()
29 : AliAnalysisTaskSE(), fESD(0), fOutputList(0), fHistCentrality(0), fHistJetArea(0), fHistJetMaxPt(0), fHistJetZ(0), fHistJetNEF(0), fHistJetPtvsCent(0), fHistJetPtM3vsCent(0), fHistLeadingJetPtvsCent(0), fHistLeadingJetPtM3vsCent(0), fHistJetAreavsCent(0), fHistJetMaxPtvsCent(0), fHistJetZvsCent(0), fHistJetNEFvsCent(0), fHistNjetvsCent(0), fHistJetPtvsNtrack(0), fHistJetAreavsNtrack(0), fHistJetMaxPtvsNtrack(0), fHistJetZvsNtrack(0), fHistJetNEFvsNtrack(0), fHistNjetvsNtrack(0), fHistDeltaRho12vsCent(0), fHistDeltaRho13vsCent(0), fHistDeltaRho23vsCent(0), fHistDeltaJetPt12vsCent(0), fHistDeltaJetPt13vsCent(0), fHistDeltaJetPt23vsCent(0), fHistRho1vsCent(0), fHistRho2vsCent(0), fHistRho3vsCent(0),
30 fTracksName("tracks"),
32 fClustersName("clusters"),
33 fRhos1Name("fArrRhos"),
34 fRhos2Name("fArrRhos"),
35 fRhos3Name("fArrRhos"),
36 phimin(-10), phimax(10),
37 etamin(-0.9), etamax(0.9),
42 for (Int_t i = 0;i<6;++i){
44 fHistAreavsRawPt[i] = 0;
45 for (Int_t j = 0;j<4;++j){
46 fHistNEFvsPt[i][j] = 0;
48 fHistZchvsPt[i][j] = 0;
49 fHistZemvsPt[i][j] = 0;
50 fHistAreavsPt[i][j] = 0;
52 fHistNconsvsPt[i][j] = 0;
53 fHistJetPt5[i][j] = 0;
54 fHistJetPt6[i][j] = 0;
55 fHistJetPt7[i][j] = 0;
56 fHistJetPt8[i][j] = 0;
60 DefineInput(0, TChain::Class());
61 DefineOutput(1, TList::Class());
63 //________________________________________________________________________
64 AliAnalysisTaskEmcalJetSpectra::AliAnalysisTaskEmcalJetSpectra(const char *name)
65 : AliAnalysisTaskSE(name), fESD(0), fOutputList(0), fHistCentrality(0), fHistJetArea(0), fHistJetMaxPt(0), fHistJetZ(0), fHistJetNEF(0), fHistJetPtvsCent(0), fHistJetPtM3vsCent(0), fHistLeadingJetPtvsCent(0), fHistLeadingJetPtM3vsCent(0), fHistJetAreavsCent(0), fHistJetMaxPtvsCent(0), fHistJetZvsCent(0), fHistJetNEFvsCent(0), fHistNjetvsCent(0), fHistJetPtvsNtrack(0), fHistJetAreavsNtrack(0), fHistJetMaxPtvsNtrack(0), fHistJetZvsNtrack(0), fHistJetNEFvsNtrack(0), fHistNjetvsNtrack(0), fHistDeltaRho12vsCent(0), fHistDeltaRho13vsCent(0), fHistDeltaRho23vsCent(0), fHistDeltaJetPt12vsCent(0), fHistDeltaJetPt13vsCent(0), fHistDeltaJetPt23vsCent(0), fHistRho1vsCent(0), fHistRho2vsCent(0), fHistRho3vsCent(0),
66 fTracksName("tracks"),
68 fClustersName("clusters"),
69 fRhos1Name("fArrRhos"),
70 fRhos2Name("fArrRhos"),
71 fRhos3Name("fArrRhos"),
72 phimin(-10), phimax(10),
73 etamin(-0.9), etamax(0.9),
77 for (Int_t i = 0;i<6;++i){
79 fHistAreavsRawPt[i] = 0;
80 for (Int_t j = 0;j<4;++j){
81 fHistNEFvsPt[i][j] = 0;
83 fHistZchvsPt[i][j] = 0;
84 fHistZemvsPt[i][j] = 0;
85 fHistAreavsPt[i][j] = 0;
87 fHistNconsvsPt[i][j] = 0;
88 fHistJetPt5[i][j] = 0;
89 fHistJetPt6[i][j] = 0;
90 fHistJetPt7[i][j] = 0;
91 fHistJetPt8[i][j] = 0;
95 DefineInput(0, TChain::Class());
96 DefineOutput(1, TList::Class());
99 //________________________________________________________________________
100 void AliAnalysisTaskEmcalJetSpectra::UserCreateOutputObjects()
103 AliVEventHandler* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
105 AliError("Input handler not available!");
110 fOutputList = new TList();
111 fOutputList->SetOwner();
113 fHistCentrality = new TH1F("Centrality", "Centrality", 101, -1, 100);
114 fHistJetPtvsCent = new TH2F("JetPtvsCent", "JetPtvsCent", 101, -1, 100, 200, 0, 500);
115 fHistJetPtM3vsCent = new TH2F("JetPtM3vsCent", "JetPtM3vsCent", 101, -1, 100, 200, -250, 250);
116 fHistLeadingJetPtvsCent = new TH2F("LeadingJetPtvsCent", "LeadingJetPtvsCent", 101, -1, 100, 200, 0, 500);
117 fHistLeadingJetPtM3vsCent = new TH2F("LeadingJetPtM3vsCent", "LeadingJetPtM3vsCent", 101, -1, 100, 200, -250, 250);
118 fHistJetAreavsCent = new TH2F("JetAreavsCent", "JetAreavsCent", 101, -1, 100, 100, 0, 1.0);
119 fHistJetMaxPtvsCent = new TH2F("JetMaxPtvsCent", "JetMaxPtvsCent", 101, -1, 100, 200, 0, 100);
120 fHistJetZvsCent = new TH2F("JetZvsCent", "JetZvsCent", 101, -1, 100, 100, 0, 1.0);
121 fHistJetNEFvsCent = new TH2F("JetNEFvsCent", "JetNEFvsCent", 101, -1, 100, 100, 0, 1.0);
122 fHistNjetvsCent = new TH2F("NjetvsCent", "NjetvsCent", 101, -1, 100, 100, 0, 100);
124 fHistJetPtvsNtrack = new TH2F("JetPtvsNtrack", "JetPtvsNtrack", 500, 0, 2500, 200, 0, 500);
125 fHistJetAreavsNtrack = new TH2F("JetAreavsNtrack", "JetAreavsNtrack", 500, 0, 2500, 100, 0, 1.0);
126 fHistJetMaxPtvsNtrack = new TH2F("JetMaxPtvsNtrack", "JetMaxPtvsNtrack", 500, 0, 2500, 200, 0, 100);
127 fHistJetZvsNtrack = new TH2F("JetZvsNtrack", "JetZvsNtrack", 500, 0, 2500, 100, 0, 1.0);
128 fHistJetNEFvsNtrack = new TH2F("JetNEFvsNtrack", "JetNEFvsNtrack", 500, 0, 2500, 100, 0, 1.0);
129 fHistNjetvsNtrack = new TH2F("NjetvsNtrack", "rNjetvsNtrack", 500, 0, 2500, 100, 0, 100);
131 fHistJetArea = new TH1F("JetArea", "Jet Area", 100, 0.0, 1.0);
132 fHistJetMaxPt = new TH1F("JetMaxPt", "Jet pt max track", 100, 0.0, 100.0);
133 fHistJetZ = new TH1F("JetZ", "Jet Z", 100, 0.0, 1.0);
134 fHistJetNEF = new TH1F("JetNEF", "Jet NEF", 100, 0.0, 1.0);
135 fHistDeltaRho12vsCent = new TH2F("DeltaRho12vsCent", "DeltaRho12vsCent", 100, 0.0, 100.0, 500,-250,250);
136 fHistDeltaRho13vsCent = new TH2F("DeltaRho13vsCent", "DeltaRho13vsCent", 100, 0.0, 100.0, 500,-250,250);
137 fHistDeltaRho23vsCent = new TH2F("DeltaRho23vsCent", "DeltaRho23vsCent", 100, 0.0, 100.0, 500,-250,250);
138 fHistDeltaJetPt12vsCent = new TH2F("DeltaJetPt12vsCent", "DeltaJetPt12vsCent", 100, 0.0, 100.0, 500,-250,250);
139 fHistDeltaJetPt13vsCent = new TH2F("DeltaJetPt13vsCent", "DeltaJetPt13vsCent", 100, 0.0, 100.0, 500,-250,250);
140 fHistDeltaJetPt23vsCent = new TH2F("DeltaJetPt23vsCent", "DeltaJetPt23vsCent", 100, 0.0, 100.0, 500,-250,250);
141 fHistRho1vsCent = new TH2F("Rho1vsCent", "Rho1vsCent", 100, 0.0, 100.0, 500, 0, 500);
142 fHistRho2vsCent = new TH2F("Rho2vsCent", "Rho2vsCent", 100, 0.0, 100.0, 500, 0, 500);
143 fHistRho3vsCent = new TH2F("Rho3vsCent", "Rho3vsCent", 100, 0.0, 100.0, 500, 0, 500);
146 for (Int_t i = 0;i<6;++i){
147 TString name00(Form("fHistRawJetPt_%i",i));
148 fHistRawJetPt[i] = new TH1F(name00,name00,250,0,500);
149 fOutputList->Add(fHistRawJetPt[i]);
150 TString name01(Form("fHistAreavsRawPt_%i",i));
151 fHistAreavsRawPt[i] = new TH2F(name01,name01,250,0,500,100,0,1);
152 fOutputList->Add(fHistAreavsRawPt[i]);
153 for (Int_t j = 0;j<4;j++){
154 TString name0(Form("fHistNEFvsPt_%i_%i",i,j));
155 fHistNEFvsPt[i][j] = new TH2F(name0,name0,250,-250,250,200,0,2);
156 fOutputList->Add(fHistNEFvsPt[i][j]);
157 TString name1(Form("fHistZvsPt_%i_%i",i,j));
158 fHistZvsPt[i][j] = new TH2F(name1,name1,250,-250,250,200,0,2);
159 fOutputList->Add(fHistZvsPt[i][j]);
160 TString name2(Form("fHistZchvsPt_%i_%i",i,j));
161 fHistZchvsPt[i][j] = new TH2F(name2,name2,250,-250,250,200,0,2);
162 fOutputList->Add(fHistZchvsPt[i][j]);
163 TString name3(Form("fHistZemvsPt_%i_%i",i,j));
164 fHistZemvsPt[i][j] = new TH2F(name3,name3,250,-250,250,200,0,2);
165 fOutputList->Add(fHistZemvsPt[i][j]);
166 TString name4(Form("fHistAreavsPt_%i_%i",i,j));
167 fHistAreavsPt[i][j] = new TH2F(name4,name4,250,-250,250,100,0,1);
168 fOutputList->Add(fHistAreavsPt[i][j]);
169 TString name5(Form("fHistJetPt_%i_%i",i,j));
170 fHistJetPt[i][j] = new TH1F(name5,name5,250,-250,250);
171 fOutputList->Add(fHistJetPt[i][j]);
172 TString name6(Form("fHistNconsvsPt_%i_%i",i,j));
173 fHistNconsvsPt[i][j] = new TH2F(name6,name6,250,-250,250,500,0,500);
174 fOutputList->Add(fHistNconsvsPt[i][j]);
175 TString name7(Form("fHistJetPt5_%i_%i",i,j));
176 fHistJetPt5[i][j] = new TH1F(name7,name7,250,-250,250);
177 fOutputList->Add(fHistJetPt5[i][j]);
178 TString name8(Form("fHistJetPt6_%i_%i",i,j));
179 fHistJetPt6[i][j] = new TH1F(name8,name8,250,-250,250);
180 fOutputList->Add(fHistJetPt6[i][j]);
181 TString name9(Form("fHistJetPt7_%i_%i",i,j));
182 fHistJetPt7[i][j] = new TH1F(name9,name9,250,-250,250);
183 fOutputList->Add(fHistJetPt7[i][j]);
184 TString name10(Form("fHistJetPt8_%i_%i",i,j));
185 fHistJetPt8[i][j] = new TH1F(name10,name10,250,-250,250);
186 fOutputList->Add(fHistJetPt8[i][j]);
189 fOutputList->Add(fHistCentrality);
190 fOutputList->Add(fHistJetArea);
191 fOutputList->Add(fHistJetMaxPt);
192 fOutputList->Add(fHistJetZ);
193 fOutputList->Add(fHistJetNEF);
194 fOutputList->Add(fHistJetPtvsCent);
195 fOutputList->Add(fHistJetPtM3vsCent);
196 fOutputList->Add(fHistLeadingJetPtvsCent);
197 fOutputList->Add(fHistLeadingJetPtM3vsCent);
198 fOutputList->Add(fHistJetAreavsCent);
199 fOutputList->Add(fHistJetMaxPtvsCent);
200 fOutputList->Add(fHistJetZvsCent);
201 fOutputList->Add(fHistJetNEFvsCent);
202 fOutputList->Add(fHistNjetvsCent);
204 fOutputList->Add(fHistJetPtvsNtrack);
205 fOutputList->Add(fHistJetAreavsNtrack);
206 fOutputList->Add(fHistJetMaxPtvsNtrack);
207 fOutputList->Add(fHistJetZvsNtrack);
208 fOutputList->Add(fHistJetNEFvsNtrack);
209 fOutputList->Add(fHistNjetvsNtrack);
210 fOutputList->Add(fHistDeltaRho12vsCent);
211 fOutputList->Add(fHistDeltaRho13vsCent);
212 fOutputList->Add(fHistDeltaRho23vsCent);
213 fOutputList->Add(fHistDeltaJetPt12vsCent);
214 fOutputList->Add(fHistDeltaJetPt13vsCent);
215 fOutputList->Add(fHistDeltaJetPt23vsCent);
216 fOutputList->Add(fHistRho1vsCent);
217 fOutputList->Add(fHistRho2vsCent);
218 fOutputList->Add(fHistRho3vsCent);
220 PostData(1, fOutputList);
223 //________________________________________________________________________
225 Int_t AliAnalysisTaskEmcalJetSpectra::GetCentBin(Double_t cent) const {
227 if (cent>=0 && cent<10)
229 else if (cent>=10 && cent<20)
231 else if (cent>=20 && cent<30)
233 else if (cent>=30 && cent<40)
235 else if (cent>=40 && cent<50)
237 else if (cent>=50 && cent<90)
241 //________________________________________________________________________
242 void AliAnalysisTaskEmcalJetSpectra::UserExec(Option_t *)
244 // Main loop, called for each event.
247 Bool_t esdMode = kTRUE;
248 if (dynamic_cast<AliAODEvent*>(InputEvent()))
251 // optimization in case autobranch loading is off
252 AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
253 if (fTracksName == "Tracks")
254 am->LoadBranch("Tracks");
258 TList *l = InputEvent()->GetList();
259 AliCentrality *centrality = InputEvent()->GetCentrality() ;
261 fCent = centrality->GetCentralityPercentile("V0M");
263 fCent=99; // probably pp data
265 AliError(Form("Centrality negative: %f", fCent));
269 Int_t centbin = GetCentBin(fCent);
271 //don't want to analyze events 90-100
275 TClonesArray *jets = 0;
276 TClonesArray *tracks = 0;
278 tracks = dynamic_cast<TClonesArray*>(l->FindObject(fTracksName));
280 AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() ));
284 jets = dynamic_cast<TClonesArray*>(l->FindObject(fJetsName));
286 AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() ));
291 TParameter<Double_t> *Rho1Param = dynamic_cast<TParameter<Double_t>*>(InputEvent()->FindListObject(fRhos1Name));
293 rho1 = Rho1Param->GetVal();
296 TParameter<Double_t> *Rho2Param = dynamic_cast<TParameter<Double_t>*>(InputEvent()->FindListObject(fRhos2Name));
298 rho2 = Rho2Param->GetVal();
301 TParameter<Double_t> *Rho3Param = dynamic_cast<TParameter<Double_t>*>(InputEvent()->FindListObject(fRhos3Name));
303 rho3 = Rho3Param->GetVal();
306 fHistRho1vsCent->Fill(fCent,rho1);
308 fHistRho2vsCent->Fill(fCent,rho2);
310 fHistRho3vsCent->Fill(fCent,rho3);
312 if (( rho1>0 ) && ( rho2>0 ))
313 fHistDeltaRho12vsCent->Fill(fCent,rho1-rho2);
314 if (( rho1>0 ) && ( rho3>0 ))
315 fHistDeltaRho13vsCent->Fill(fCent,rho1-rho3);
316 if (( rho2>0 ) && ( rho3>0 ))
317 fHistDeltaRho23vsCent->Fill(fCent,rho2-rho3);
318 fHistCentrality->Fill(fCent);
320 const Int_t Ntracks = tracks->GetEntries();
321 Double_t lJetPt = -500;
322 Double_t lJetPtun = -500;
324 const Int_t Njets = jets->GetEntries();
328 for (Int_t iJets = 0; iJets < Njets; ++iJets) {
330 AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(iJets));
333 if ((jet->Phi()<phimin)||(jet->Phi()>phimax))
335 if ((jet->Eta()<etamin)||(jet->Eta()>etamax))
337 fHistAreavsRawPt[centbin]->Fill(jet->Pt(),jet->Area());
338 if (jet->Area()<areacut)
340 //prevents 0 area jets from sneaking by when area cut == 0
343 fHistRawJetPt[centbin]->Fill(jet->Pt());
346 fHistJetArea->Fill(jet->Area());
347 fHistJetMaxPt->Fill(jet->MaxTrackPt());
353 jetPt1 = jet->Pt()-jet->Area()*rho1;
355 jetPt2 = jet->Pt()-jet->Area()*rho2;
357 jetPt3 = jet->Pt()-jet->Area()*rho3;
358 if (( rho1>0 ) && ( rho2>0 ))
359 fHistDeltaJetPt12vsCent->Fill(fCent,jetPt1-jetPt2);
360 if (( rho1>0 ) && ( rho3>0 ))
361 fHistDeltaJetPt13vsCent->Fill(fCent,jetPt1-jetPt3);
362 if (( rho2>0 ) && ( rho3>0 ))
363 fHistDeltaJetPt23vsCent->Fill(fCent,jetPt2-jetPt3);
365 if (lJetPt < jetPt2){
367 lJetPtun = jet->Pt();
370 if (jet->MaxTrackPt()>jet->MaxClusterPt()){
371 Z = jet->MaxTrackPt();
372 fHistJetMaxPtvsCent->Fill(fCent,jet->MaxTrackPt());
373 fHistJetMaxPtvsNtrack->Fill(Ntracks,jet->MaxTrackPt());
376 Z = jet->MaxClusterPt();
377 fHistJetMaxPtvsCent->Fill(fCent,jet->MaxClusterPt());
378 fHistJetMaxPtvsNtrack->Fill(Ntracks,jet->MaxClusterPt());
381 fHistJetZ->Fill(jet->MaxTrackPt()/jetPt2);
383 fHistJetNEF->Fill(jet->NEF());
384 fHistJetPtvsCent->Fill(fCent,jet->Pt());
385 fHistJetPtM3vsCent->Fill(fCent,jetPt2);
386 fHistJetAreavsCent->Fill(fCent,jet->Area());
387 fHistJetZvsCent->Fill(fCent,Z);
388 fHistJetZvsNtrack->Fill(Ntracks,Z);
390 fHistNEFvsPt[centbin][0]->Fill(jet->Pt(),jet->NEF());
391 fHistZvsPt[centbin][0]->Fill(jet->Pt(),Z/jet->Pt());
392 fHistZchvsPt[centbin][0]->Fill(jet->Pt(),jet->MaxTrackPt()/jet->Pt());
393 fHistZemvsPt[centbin][0]->Fill(jet->Pt(),jet->MaxClusterPt()/jet->Pt());
394 fHistAreavsPt[centbin][0]->Fill(jet->Pt(),jet->Area());
395 fHistJetPt[centbin][0]->Fill(jet->Pt());
396 fHistNconsvsPt[centbin][0]->Fill(jet->Pt(),jet->N());
397 if ((jet->MaxTrackPt()>5) || (jet->MaxClusterPt()>5))
398 fHistJetPt5[centbin][0]->Fill(jet->Pt());
399 if ((jet->MaxTrackPt()>6) || (jet->MaxClusterPt()>6))
400 fHistJetPt6[centbin][0]->Fill(jet->Pt());
401 if ((jet->MaxTrackPt()>7) || (jet->MaxClusterPt()>7))
402 fHistJetPt7[centbin][0]->Fill(jet->Pt());
403 if ((jet->MaxTrackPt()>8) || (jet->MaxClusterPt()>8))
404 fHistJetPt8[centbin][0]->Fill(jet->Pt());
406 fHistNEFvsPt[centbin][1]->Fill(jetPt1,jet->NEF());
407 fHistZvsPt[centbin][1]->Fill(jetPt1,Z/jetPt1);
408 fHistZchvsPt[centbin][1]->Fill(jetPt1,jet->MaxTrackPt()/jetPt1);
409 fHistZemvsPt[centbin][1]->Fill(jetPt1,jet->MaxClusterPt()/jetPt1);
410 fHistAreavsPt[centbin][1]->Fill(jetPt1,jet->Area());
411 fHistJetPt[centbin][1]->Fill(jetPt1);
412 fHistNconsvsPt[centbin][1]->Fill(jetPt1,jet->N());
413 if ((jet->MaxTrackPt()>5) || (jet->MaxClusterPt()>5))
414 fHistJetPt5[centbin][1]->Fill(jetPt1);
415 if ((jet->MaxTrackPt()>6) || (jet->MaxClusterPt()>6))
416 fHistJetPt6[centbin][1]->Fill(jetPt1);
417 if ((jet->MaxTrackPt()>7) || (jet->MaxClusterPt()>7))
418 fHistJetPt7[centbin][1]->Fill(jetPt1);
419 if ((jet->MaxTrackPt()>8) || (jet->MaxClusterPt()>8))
420 fHistJetPt8[centbin][1]->Fill(jetPt1);
422 fHistNEFvsPt[centbin][2]->Fill(jetPt2,jet->NEF());
423 fHistZvsPt[centbin][2]->Fill(jetPt2,Z/jetPt2);
424 fHistZchvsPt[centbin][2]->Fill(jetPt2,jet->MaxTrackPt()/jetPt2);
425 fHistZemvsPt[centbin][2]->Fill(jetPt2,jet->MaxClusterPt()/jetPt2);
426 fHistAreavsPt[centbin][2]->Fill(jetPt2,jet->Area());
427 fHistJetPt[centbin][2]->Fill(jetPt2);
428 fHistNconsvsPt[centbin][2]->Fill(jetPt2,jet->N());
429 if ((jet->MaxTrackPt()>5) || (jet->MaxClusterPt()>5))
430 fHistJetPt5[centbin][2]->Fill(jetPt2);
431 if ((jet->MaxTrackPt()>6) || (jet->MaxClusterPt()>6))
432 fHistJetPt6[centbin][2]->Fill(jetPt2);
433 if ((jet->MaxTrackPt()>7) || (jet->MaxClusterPt()>7))
434 fHistJetPt7[centbin][2]->Fill(jetPt2);
435 if ((jet->MaxTrackPt()>8) || (jet->MaxClusterPt()>8))
436 fHistJetPt8[centbin][2]->Fill(jetPt2);
438 fHistNEFvsPt[centbin][3]->Fill(jetPt3,jet->NEF());
439 fHistZvsPt[centbin][3]->Fill(jetPt3,Z/jetPt3);
440 fHistZchvsPt[centbin][3]->Fill(jetPt3,jet->MaxTrackPt()/jetPt3);
441 fHistZemvsPt[centbin][3]->Fill(jetPt3,jet->MaxClusterPt()/jetPt3);
442 fHistAreavsPt[centbin][3]->Fill(jetPt3,jet->Area());
443 fHistJetPt[centbin][3]->Fill(jetPt3);
444 fHistNconsvsPt[centbin][3]->Fill(jetPt3,jet->N());
445 if ((jet->MaxTrackPt()>5) || (jet->MaxClusterPt()>5))
446 fHistJetPt5[centbin][3]->Fill(jetPt3);
447 if ((jet->MaxTrackPt()>6) || (jet->MaxClusterPt()>6))
448 fHistJetPt6[centbin][3]->Fill(jetPt3);
449 if ((jet->MaxTrackPt()>7) || (jet->MaxClusterPt()>7))
450 fHistJetPt7[centbin][3]->Fill(jetPt3);
451 if ((jet->MaxTrackPt()>8) || (jet->MaxClusterPt()>8))
452 fHistJetPt8[centbin][3]->Fill(jetPt3);
454 fHistJetNEFvsCent->Fill(fCent,jet->NEF());
455 fHistJetPtvsNtrack->Fill(Ntracks,jetPt2);
456 fHistJetAreavsNtrack->Fill(Ntracks,jet->Area());
457 fHistJetNEFvsNtrack->Fill(Ntracks,jet->NEF());
460 fHistLeadingJetPtvsCent->Fill(fCent,lJetPtun);
461 fHistLeadingJetPtM3vsCent->Fill(fCent,lJetPt);
462 fHistNjetvsCent->Fill(fCent,NjetAcc);
463 fHistNjetvsNtrack->Fill(Ntracks,NjetAcc);
465 PostData(1, fOutputList);
468 //________________________________________________________________________
469 void AliAnalysisTaskEmcalJetSpectra::Terminate(Option_t *)