]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetHMEC.cxx
Charged jets (pPb): Improved trackcut analysis
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskEmcalJetHMEC.cxx
CommitLineData
7f71d351 1// $Id$
5b38fc6d 2
1af24178 3//////////
4//Measure Jet-hadron correlations
5//Does event Mixing using AliEventPoolManager
6/////////
7
5b38fc6d 8#include "AliAnalysisTaskEmcalJetHMEC.h"
9
1101e468 10#include "TChain.h"
11#include "TTree.h"
12#include "TList.h"
13#include "TH1F.h"
14#include "TH2F.h"
15#include "THnSparse.h"
16#include "TCanvas.h"
17#include <TClonesArray.h>
18#include <TParticle.h>
19#include "AliVTrack.h"
20#include "TParameter.h"
21
22#include "AliAODEvent.h"
23#include "AliAnalysisTask.h"
24#include "AliAnalysisManager.h"
25
26#include "AliESDEvent.h"
27#include "AliESDInputHandler.h"
28#include "AliESDCaloCluster.h"
29#include "AliESDVertex.h"
30#include "AliCentrality.h"
31#include "AliAODJet.h"
32#include "AliEmcalJet.h"
33#include "AliESDtrackCuts.h"
1101e468 34
5b38fc6d 35#include "TVector3.h"
36#include "AliPicoTrack.h"
55f64c2f 37#include "AliEventPoolManager.h"
5b38fc6d 38
1101e468 39ClassImp(AliAnalysisTaskEmcalJetHMEC)
40
41//________________________________________________________________________
42AliAnalysisTaskEmcalJetHMEC::AliAnalysisTaskEmcalJetHMEC() :
f569a5d2 43 AliAnalysisTaskEmcalJet("HMEC",kFALSE),
44 fTracksName(""),
45 fJetsName(""),
1101e468 46 fPhimin(-10),
47 fPhimax(10),
48 fEtamin(-0.9),
49 fEtamax(0.9),
50 fAreacut(0.0),
965c985f 51 fTrkBias(5),
52 fClusBias(5),
22d3f9c2 53 fTrkEta(0.9),
fab37389 54 fDoEventMixing(0),
55 fMixingTracks(50000),
f569a5d2 56 fESD(0),
57 fAOD(0),
55f64c2f 58 fPoolMgr(0x0),
1101e468 59 fHistTrackPt(0),
60 fHistCentrality(0),
61 fHistJetEtaPhi(0),
fab37389 62 fHistJetHEtaPhi(0),
22d3f9c2 63 fhnMixedEvents(0x0),
64 fhnJH(0x0)
1101e468 65{
66 // Default Constructor
67
22d3f9c2 68 for(Int_t ipta=0; ipta<7; ipta++){
69 fHistTrackEtaPhi[ipta]=0;
70 }
71
1101e468 72 for(Int_t icent = 0; icent<6; ++icent){
73 fHistJetPt[icent]=0;
81941f29 74 fHistJetPtBias[icent]=0;
41f87d59 75 fHistLeadJetPt[icent]=0;
76 fHistLeadJetPtBias[icent]=0;
81941f29 77 fHistJetPtTT[icent]=0;
55f64c2f 78 for(Int_t iptjet = 0; iptjet<5; ++iptjet){
1101e468 79 for(Int_t ieta = 0; ieta<3; ++ieta){
f569a5d2 80 fHistJetH[icent][iptjet][ieta]=0;
81 fHistJetHBias[icent][iptjet][ieta]=0;
82 fHistJetHTT[icent][iptjet][ieta]=0;
1101e468 83 }
84 }
85 }
86
87}
88//________________________________________________________________________
89AliAnalysisTaskEmcalJetHMEC::AliAnalysisTaskEmcalJetHMEC(const char *name) :
f569a5d2 90 AliAnalysisTaskEmcalJet(name,kTRUE),
91 fTracksName(""),
92 fJetsName(""),
1101e468 93 fPhimin(-10),
94 fPhimax(10),
95 fEtamin(-0.9),
96 fEtamax(0.9),
97 fAreacut(0.0),
965c985f 98 fTrkBias(5),
99 fClusBias(5),
22d3f9c2 100 fTrkEta(0.9),
55f64c2f 101 fDoEventMixing(0),
102 fMixingTracks(50000),
f569a5d2 103 fESD(0),
104 fAOD(0),
55f64c2f 105 fPoolMgr(0x0),
1101e468 106 fHistTrackPt(0),
107 fHistCentrality(0),
108 fHistJetEtaPhi(0),
55f64c2f 109 fHistJetHEtaPhi(0),
22d3f9c2 110 fhnMixedEvents(0x0),
111 fhnJH(0x0)
1101e468 112{
113 // Constructor
22d3f9c2 114 for(Int_t ipta=0; ipta<7; ipta++){
115 fHistTrackEtaPhi[ipta]=0;
116 }
1101e468 117 for(Int_t icent = 0; icent<6; ++icent){
118 fHistJetPt[icent]=0;
119 fHistJetPtBias[icent]=0;
41f87d59 120 fHistLeadJetPt[icent]=0;
121 fHistLeadJetPtBias[icent]=0;
81941f29 122 fHistJetPtTT[icent]=0;
123 for(Int_t iptjet = 0; iptjet<5; ++iptjet){
1101e468 124 for(Int_t ieta = 0; ieta<3; ++ieta){
f569a5d2 125 fHistJetH[icent][iptjet][ieta]=0;
126 fHistJetHBias[icent][iptjet][ieta]=0;
127 fHistJetHTT[icent][iptjet][ieta]=0;
1101e468 128 }
129 }
130 }
131
1101e468 132}
133
134//________________________________________________________________________
f569a5d2 135void AliAnalysisTaskEmcalJetHMEC::UserCreateOutputObjects() {
136 // Called once
137 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
1101e468 138 OpenFile(1);
1101e468 139
55f64c2f 140 // Create histograms
1101e468 141 fHistTrackPt = new TH1F("fHistTrackPt", "P_{T} distribution", 1000, 0.0, 100.0);
1101e468 142 fHistCentrality = new TH1F("fHistCentrality","centrality",100,0,100);
9e4a0350 143 fHistJetEtaPhi = new TH2F("fHistJetEtaPhi","Jet eta-phi",900,-1.8,1.8,720,-3.2,3.2);
144 fHistJetHEtaPhi = new TH2F("fHistJetHEtaPhi","Jet-Hadron deta-dphi",900,-1.8,1.8,720,-1.6,4.8);
1101e468 145
daf7dbf5 146 TString name;
1101e468 147
22d3f9c2 148 for(Int_t ipta=0; ipta<7; ++ipta){
daf7dbf5 149 name = Form("fHistTrackEtaPhi_%i", ipta);
9e4a0350 150 fHistTrackEtaPhi[ipta] = new TH2F(name,name,400,-1,1,720,0.0,2.0*TMath::Pi());
f569a5d2 151 fOutput->Add(fHistTrackEtaPhi[ipta]);
22d3f9c2 152 }
1101e468 153
154 for(Int_t icent = 0; icent<6; ++icent){
daf7dbf5 155 name = Form("fHistJetPt_%i",icent);
1101e468 156 fHistJetPt[icent] = new TH1F(name,name,200,0,200);
f569a5d2 157 fOutput->Add(fHistJetPt[icent]);
1101e468 158
daf7dbf5 159 name = Form("fHistJetPtBias_%i",icent);
1101e468 160 fHistJetPtBias[icent] = new TH1F(name,name,200,0,200);
f569a5d2 161 fOutput->Add(fHistJetPtBias[icent]);
1101e468 162
daf7dbf5 163 name = Form("fHistLeadJetPt_%i",icent);
41f87d59 164 fHistLeadJetPt[icent] = new TH1F(name,name,200,0,200);
f569a5d2 165 fOutput->Add(fHistLeadJetPt[icent]);
41f87d59 166
daf7dbf5 167 name = Form("fHistLeadJetPtBias_%i",icent);
41f87d59 168 fHistLeadJetPtBias[icent] = new TH1F(name,name,200,0,200);
f569a5d2 169 fOutput->Add(fHistLeadJetPtBias[icent]);
41f87d59 170
daf7dbf5 171 name = Form("fHistJetPtTT_%i",icent);
81941f29 172 fHistJetPtTT[icent] = new TH1F(name,name,200,0,200);
f569a5d2 173 fOutput->Add(fHistJetPtTT[icent]);
81941f29 174
175 for(Int_t iptjet = 0; iptjet<5; ++iptjet){
1101e468 176 for(Int_t ieta = 0; ieta<3; ++ieta){
f569a5d2 177 name = Form("fHistJetH_%i_%i_%i",icent,iptjet,ieta);
178 fHistJetH[icent][iptjet][ieta]=new TH2F(name,name,72,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30);
179 fOutput->Add(fHistJetH[icent][iptjet][ieta]);
1101e468 180
f569a5d2 181 name = Form("fHistJetHBias_%i_%i_%i",icent,iptjet,ieta);
182 fHistJetHBias[icent][iptjet][ieta]=new TH2F(name,name,72,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30);
183 fOutput->Add(fHistJetHBias[icent][iptjet][ieta]);
1101e468 184
f569a5d2 185 name = Form("fHistJetHTT_%i_%i_%i",icent,iptjet,ieta);
186 fHistJetHTT[icent][iptjet][ieta]=new TH2F(name,name,72,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30);
187 fOutput->Add(fHistJetHTT[icent][iptjet][ieta]);
81941f29 188
1101e468 189 }
190 }
191 }
192
22d3f9c2 193 UInt_t cifras = 0; // bit coded, see GetDimParams() below
c1076e58 194 cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7 | 1<<8;
22d3f9c2 195 fhnJH = NewTHnSparseF("fhnJH", cifras);
22d3f9c2 196 fhnJH->Sumw2();
f569a5d2 197 fOutput->Add(fhnJH);
22d3f9c2 198
55f64c2f 199 if(fDoEventMixing){
c1076e58 200 cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7 | 1<<8;
201 fhnMixedEvents = NewTHnSparseF("fhnMixedEvents", cifras);
c1076e58 202 fhnMixedEvents->Sumw2();
f569a5d2 203 fOutput->Add(fhnMixedEvents);
965c985f 204 }
205
f569a5d2 206 fOutput->Add(fHistTrackPt);
207 fOutput->Add(fHistCentrality);
208 fOutput->Add(fHistJetEtaPhi);
209 fOutput->Add(fHistJetHEtaPhi);
1101e468 210
f569a5d2 211 PostData(1, fOutput);
55f64c2f 212
213 //Event Mixing
214 Int_t trackDepth = fMixingTracks;
215 Int_t poolsize = 1000; // Maximum number of events, ignored in the present implemented of AliEventPoolManager
216
dc64e292 217 Int_t nZvtxBins = 5+1+5;
55f64c2f 218 // bins for second buffer are shifted by 100 cm
dc64e292 219 Double_t vertexBins[] = { -10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10, };
55f64c2f 220 Double_t* zvtxbin = vertexBins;
221
222 Int_t nCentralityBins = 100;
223 Double_t centralityBins[nCentralityBins];
224 for(Int_t ic=0; ic<nCentralityBins; ic++){
225 centralityBins[ic]=1.0*ic;
226 }
55f64c2f 227
228 fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBins, centralityBins, nZvtxBins, zvtxbin);
229
1101e468 230}
231
232//________________________________________________________________________
233
234Double_t AliAnalysisTaskEmcalJetHMEC:: RelativePhi(Double_t mphi,Double_t vphi) {
235
1101e468 236 if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
237 else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
238 if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
239 else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
240 Double_t dphi = vphi-mphi;
241 if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
242 else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
243
244 return dphi;//dphi in [-Pi, Pi]
245}
246
1101e468 247//________________________________________________________________________
f569a5d2 248Int_t AliAnalysisTaskEmcalJetHMEC::GetCentBin(Double_t cent) const {
1101e468 249 // Get centrality bin.
250
251 Int_t centbin = -1;
f569a5d2 252 if (cent>=0 && cent<10) centbin = 0;
253 else if (cent>=10 && cent<20) centbin = 1;
254 else if (cent>=20 && cent<30) centbin = 2;
255 else if (cent>=30 && cent<40) centbin = 3;
256 else if (cent>=40 && cent<50) centbin = 4;
257 else if (cent>=50 && cent<90) centbin = 5;
1101e468 258 return centbin;
259}
260
1101e468 261//________________________________________________________________________
f569a5d2 262Int_t AliAnalysisTaskEmcalJetHMEC::GetEtaBin(Double_t eta) const {
1101e468 263 // Get eta bin for histos.
264
265 Int_t etabin = -1;
f569a5d2 266 if (TMath::Abs(eta)<=0.4) etabin = 0;
267 else if (TMath::Abs(eta)>0.4 && TMath::Abs(eta)<0.8) etabin = 1;
268 else if (TMath::Abs(eta)>=0.8) etabin = 2;
1101e468 269 return etabin;
270}
271//________________________________________________________________________
81941f29 272Int_t AliAnalysisTaskEmcalJetHMEC::GetpTjetBin(Double_t pt) const
1101e468 273{
81941f29 274 // Get jet pt bin for histos.
275
276 Int_t ptbin = -1;
f569a5d2 277 if (pt>=15 && pt<20) ptbin = 0;
278 else if (pt>=20 && pt<25) ptbin = 1;
279 else if (pt>=25 && pt<30) ptbin = 2;
280 else if (pt>=30 && pt<60) ptbin = 3;
281 else if (pt>=60) ptbin = 4;
81941f29 282
283 return ptbin;
1101e468 284}
285
f569a5d2 286//________________________________________________________________________
287void AliAnalysisTaskEmcalJetHMEC::ExecOnce() {
288 AliAnalysisTaskEmcalJet::ExecOnce();
289
290}
1101e468 291
292//________________________________________________________________________
f569a5d2 293Bool_t AliAnalysisTaskEmcalJetHMEC::Run() {
294 // Main loop called for each event
295 if(!fTracks){
296 AliError(Form("No fTracks object!!\n"));
297 return kTRUE;
298 }
299 if(!fJets){
300 AliError(Form("No fJets object!!\n"));
301 return kTRUE;
302 }
1101e468 303
f569a5d2 304 // what kind of event do we have: AOD or ESD?
305 Bool_t esdMode = kTRUE;
306 if (dynamic_cast<AliAODEvent*>(InputEvent())) esdMode = kFALSE;
1101e468 307
f569a5d2 308 // if we have ESD event, set up ESD object
309 if(esdMode){
310 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
311 if (!fESD) {
312 AliError(Form("ERROR: fESD not available\n"));
313 return kTRUE;
314 }
315 }
1101e468 316
f569a5d2 317 // if we have AOD event, set up AOD object
318 if(!esdMode){
319 fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
320 if(!fAOD) {
321 AliError(Form("ERROR: fAOD not available\n"));
322 return kTRUE;
323 }
324 }
1101e468 325
f569a5d2 326 TList *list = InputEvent()->GetList();
327 if(!list) {
328 AliError(Form("ERROR: list not attached\n"));
329 return kTRUE;
1101e468 330 }
1101e468 331
f569a5d2 332 // get centrality
333 if (fCent<0) {
334 AliError(Form("Centrality negative: %f", fCent));
335 return kTRUE;
1101e468 336 }
337
f569a5d2 338 Int_t centbin = GetCentBin(fCent);
339 if(centbin<0) return kTRUE;
340
9e4a0350 341 Double_t fvertex[3]={0,0,0};
342 InputEvent()->GetPrimaryVertex()->GetXYZ(fvertex);
343 Double_t zVtx=fvertex[2];
344
f569a5d2 345 if(fabs(zVtx)>10.0) return kTRUE;
1101e468 346
f569a5d2 347 fHistCentrality->Fill(fCent);
1af24178 348
1101e468 349 TClonesArray *jets = 0;
350 TClonesArray *tracks = 0;
351
f569a5d2 352 tracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracks));
1101e468 353 if (!tracks) {
f569a5d2 354 AliError(Form("Pointer to tracks %s == 0", fTracks->GetName() ));
355 return kTRUE;
1101e468 356 }
357 const Int_t Ntracks=tracks->GetEntries();
358
f569a5d2 359 jets= dynamic_cast<TClonesArray*>(list->FindObject(fJets));
9e4a0350 360 if (!jets) {
f569a5d2 361 AliError(Form("Pointer to tracks %s == 0", fJets->GetName() ));
362 return kTRUE;
9e4a0350 363 }
1101e468 364 const Int_t Njets = jets->GetEntries();
365
81941f29 366 //Leticia's loop to find hardest track
81941f29 367 Int_t iTT=-1;
368 Double_t ptmax=-10;
369
f569a5d2 370 for (Int_t iTracks = 0; iTracks < Ntracks; iTracks++) {
371 AliVTrack* track = static_cast<AliVTrack*>(tracks->At(iTracks));
372 if (!track) {
373 printf("ERROR: Could not receive track %d\n", iTracks);
374 continue;
375 }
81941f29 376
f569a5d2 377 if(TMath::Abs(track->Eta())>0.9) continue;
378 if(track->Pt()<0.15) continue;
379 //iCount++;
380 if(track->Pt()>ptmax){
381 ptmax=track->Pt();
382 iTT=iTracks;
81941f29 383 }
f569a5d2 384 }
1101e468 385
81941f29 386 Int_t ijethi=-1;
81941f29 387 Double_t highestjetpt=0.0;
81941f29 388 Int_t passedTTcut=0;
1101e468 389
f569a5d2 390 for (Int_t ijets = 0; ijets < Njets; ijets++){
391 AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijets));
392 if (!jet) continue;
393 if(!AcceptthisJet(jet)) continue;
1101e468 394
1101e468 395 Double_t jetPt = jet->Pt();
396
59d5ae99 397 if(highestjetpt<jetPt){
f569a5d2 398 ijethi=ijets;
399 highestjetpt=jetPt;
59d5ae99 400 }
f569a5d2 401 }
59d5ae99 402
f569a5d2 403 for (Int_t ijet = 0; ijet < Njets; ijet++){
c1076e58 404 AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet));
f569a5d2 405 if (!jet) continue;
406 if(!AcceptthisJet(jet)) continue;
c1076e58 407
55f64c2f 408 Double_t jetphi = jet->Phi();
59d5ae99 409 Double_t jetPt = jet->Pt();
55f64c2f 410 Double_t jeteta=jet->Eta();
1101e468 411
41f87d59 412 Double_t leadjet=0;
413 if (ijet==ijethi) leadjet=1;
414
965c985f 415 fHistJetPt[centbin]->Fill(jet->Pt());
41f87d59 416 fHistLeadJetPt[centbin]->Fill(jet->Pt());
965c985f 417
41f87d59 418 if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){
965c985f 419 fHistJetPtBias[centbin]->Fill(jet->Pt());
41f87d59 420 fHistLeadJetPtBias[centbin]->Fill(jet->Pt());
421 }
81941f29 422
f569a5d2 423 fHistJetEtaPhi->Fill(jet->Eta(),jetphi);
1101e468 424
f569a5d2 425 if(iTT>0){
426 AliVTrack* TT = static_cast<AliVTrack*>(tracks->At(iTT));
427 if(TMath::Abs(jetphi-TT->Phi()-TMath::Pi())<0.6) passedTTcut=1;
428 else passedTTcut=0;
429 }
1101e468 430
f569a5d2 431 if(passedTTcut)
432 fHistJetPtTT[centbin]->Fill(jet->Pt());
1101e468 433
daf7dbf5 434 Int_t iptjet=-1;
435 iptjet=GetpTjetBin(jetPt);
436 if(iptjet<0) continue;
81941f29 437
f569a5d2 438 //if (highestjetpt>15) {
439 if (highestjetpt>8) {
daf7dbf5 440
f569a5d2 441 for (Int_t iTracks = 0; iTracks < Ntracks; iTracks++) {
442 AliVTrack* track = static_cast<AliVTrack*>(tracks->At(iTracks));
443 if (!track) {
444 printf("ERROR: Could not receive track %d\n", iTracks);
445 continue;
446 }
81941f29 447
f569a5d2 448 if(TMath::Abs(track->Eta())>fTrkEta) continue;
1101e468 449
f569a5d2 450 fHistTrackPt->Fill(track->Pt());
1101e468 451
f569a5d2 452 if (track->Pt()<0.15) continue;
1101e468 453
f569a5d2 454 Double_t trackphi = track->Phi();
455 if (trackphi > TMath::Pi())
456 trackphi = trackphi-2*TMath::Pi();
457
458 Double_t tracketa=track->Eta();
459 Double_t trackpt=track->Pt();
460 Double_t deta=tracketa-jeteta;
461 Int_t ieta=GetEtaBin(deta);
462 if (ieta<0) {
463 AliError(Form("Eta Bin negative: %f", deta));
464 continue;
465 }
daf7dbf5 466
f569a5d2 467 //Jet pt, track pt, dPhi,deta,fCent
468 Double_t dphijh = RelativePhi(jetphi,trackphi);
469 if (dphijh < -0.5*TMath::Pi())
470 dphijh+= 2*TMath::Pi();
471 if (dphijh > 1.5*TMath::Pi())
472 dphijh-=2.*TMath::Pi();
1101e468 473
f569a5d2 474 fHistJetH[centbin][iptjet][ieta]->Fill(dphijh,track->Pt());
475 fHistJetHEtaPhi->Fill(deta,dphijh);
22d3f9c2 476
f569a5d2 477 Double_t dR=sqrt(deta*deta+dphijh*dphijh);
1101e468 478
f569a5d2 479 if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){
480 fHistJetHBias[centbin][iptjet][ieta]->Fill(dphijh,trackpt);
22d3f9c2 481
f569a5d2 482 Double_t triggerEntries[8] = {fCent,jetPt,trackpt,dR,deta,dphijh,0.0,leadjet};
483 fhnJH->Fill(triggerEntries);
484 }
22d3f9c2 485
f569a5d2 486 if(passedTTcut)
487 fHistJetHTT[centbin][iptjet][ieta]->Fill(dphijh,trackpt);
55f64c2f 488
f569a5d2 489 } //track loop
490 }//jet pt cut
c1076e58 491 }//jet loop
55f64c2f 492
1af24178 493 //Prepare to do event mixing
494
55f64c2f 495 // create a list of reduced objects. This speeds up processing and reduces memory consumption for the event pool
496 TObjArray* tracksClone = CloneAndReduceTrackList(tracks);
497 //delete tracks;
498
55f64c2f 499 if(fDoEventMixing>0){
500
501 // event mixing
502
503 // 1. First get an event pool corresponding in mult (cent) and
504 // zvertex to the current event. Once initialized, the pool
505 // should contain nMix (reduced) events. This routine does not
506 // pre-scan the chain. The first several events of every chain
507 // will be skipped until the needed pools are filled to the
508 // specified depth. If the pool categories are not too rare, this
509 // should not be a problem. If they are rare, you could lose
510 // statistics.
511
512 // 2. Collect the whole pool's content of tracks into one TObjArray
513 // (bgTracks), which is effectively a single background super-event.
514
515 // 3. The reduced and bgTracks arrays must both be passed into
516 // FillCorrelations(). Also nMix should be passed in, so a weight
517 // of 1./nMix can be applied.
518
f569a5d2 519 AliEventPool* pool = fPoolMgr->GetEventPool(fCent, zVtx);
55f64c2f 520
daf7dbf5 521 if (!pool){
f569a5d2 522 AliFatal(Form("No pool found for centrality = %f, zVtx = %f", fCent, zVtx));
523 return kTRUE;
daf7dbf5 524 }
55f64c2f 525
1af24178 526 //check for a trigger jet
f569a5d2 527 if (pool->IsReady() || pool->NTracksInPool() > fMixingTracks / 10 || pool->GetCurrentNEvents() >= 5) {
1af24178 528
c1076e58 529 for (Int_t ijet = 0; ijet < Njets; ijet++){
c1076e58 530 Double_t leadjet=0;
531 if (ijet==ijethi) leadjet=1;
1af24178 532
c1076e58 533 AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet));
f569a5d2 534 if(!AcceptthisJet(jet)) continue;
c1076e58 535
536 Double_t jetPt = jet->Pt();
1af24178 537 Double_t jetphi = jet->Phi();
1af24178 538 Double_t jeteta=jet->Eta();
539
1af24178 540 Int_t nMix = pool->GetCurrentNEvents();
541
965c985f 542 //Fill for biased jet triggers only
543 if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){
544
545 // Fill mixed-event histos here
f569a5d2 546 for (Int_t jMix=0; jMix<nMix; jMix++) {
547 TObjArray* bgTracks = pool->GetEvent(jMix);
548 const Int_t Nbgtrks = bgTracks->GetEntries();
965c985f 549
f569a5d2 550 for(Int_t ibg=0; ibg<Nbgtrks; ibg++){
551 AliPicoTrack *part = static_cast<AliPicoTrack*>(bgTracks->At(ibg));
552 if(!part) continue;
553
554 Double_t DEta = part->Eta()-jeteta;
555 Double_t DPhi = RelativePhi(jetphi,part->Phi());
556
557 Double_t DR=TMath::Sqrt(DPhi*DPhi+DEta*DEta);
558 if(DPhi<-0.5*TMath::Pi()) DPhi+=2.*TMath::Pi();
559 if(DPhi>3./2.*TMath::Pi()) DPhi-=2.*TMath::Pi();
560 Double_t triggerEntries[8] = {fCent,jetPt,part->Pt(),DR,DEta,DPhi,0.0,leadjet};
561 fhnMixedEvents->Fill(triggerEntries,1./nMix);
562 }
1af24178 563 }
f569a5d2 564 }
965c985f 565 }
55f64c2f 566 }
55f64c2f 567
1af24178 568 //update pool if jet in event or not
569 pool->UpdatePool(tracksClone);
59d5ae99 570 }
571
f569a5d2 572 return kTRUE;
1101e468 573}
574
575//________________________________________________________________________
576void AliAnalysisTaskEmcalJetHMEC::Terminate(Option_t *)
577{
578 //just terminate
1101e468 579}
580
c1076e58 581//________________________________________________________________________
f569a5d2 582Int_t AliAnalysisTaskEmcalJetHMEC::AcceptthisJet(AliEmcalJet *jet)
c1076e58 583{
584 //applies all jet cuts except pt
585 float jetphi = jet->Phi();
586 if (jetphi>TMath::Pi())
587 jetphi = jetphi-2*TMath::Pi();
f569a5d2 588
c1076e58 589 if ((jet->Phi()<fPhimin)||(jet->Phi()>fPhimax))
590 return 0;
591 if ((jet->Eta()<fEtamin)||(jet->Eta()>fEtamax))
592 return 0;
593 if (jet->Area()<fAreacut)
594 return 0;
595 //prevents 0 area jets from sneaking by when area cut == 0
596 if (jet->Area()==0)
f569a5d2 597 return 0;
c1076e58 598 //exclude jets with extremely high pt tracks which are likely misreconstructed
599 if(jet->MaxTrackPt()>100)
600 return 0;
601
602 //passed all above cuts
603 return 1;
c1076e58 604}
605
606//________________________________________________________________________
607
f569a5d2 608THnSparse* AliAnalysisTaskEmcalJetHMEC::NewTHnSparseF(const char* name, UInt_t entries){
55f64c2f 609 // generate new THnSparseF, axes are defined in GetDimParams()
610
611 Int_t count = 0;
612 UInt_t tmp = entries;
613 while(tmp!=0){
614 count++;
615 tmp = tmp &~ -tmp; // clear lowest bit
616 }
617
618 TString hnTitle(name);
619 const Int_t dim = count;
620 Int_t nbins[dim];
621 Double_t xmin[dim];
622 Double_t xmax[dim];
623
624 Int_t i=0;
625 Int_t c=0;
626 while(c<dim && i<32){
627 if(entries&(1<<i)){
628
629 TString label("");
630 GetDimParams(i, label, nbins[c], xmin[c], xmax[c]);
631 hnTitle += Form(";%s",label.Data());
632 c++;
633 }
634
635 i++;
636 }
637 hnTitle += ";";
638
639 return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
640}
641
642void AliAnalysisTaskEmcalJetHMEC::GetDimParams(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
643{
644 // stores label and binning of axis for THnSparse
645
646 const Double_t pi = TMath::Pi();
647
648 switch(iEntry){
649
650 case 0:
651 label = "V0 centrality (%)";
55f64c2f 652 nbins = 10;
653 xmin = 0.;
654 xmax = 100.;
655 break;
656
55f64c2f 657 case 1:
658 label = "corrected jet pt";
f569a5d2 659 nbins = 20;
660 xmin = 0.;
661 xmax = 200.;
662 break;
55f64c2f 663
664 case 2:
665 label = "track pT";
f569a5d2 666 nbins = 100;
667 xmin = 0.;
668 xmax = 10;
669 break;
670
55f64c2f 671 case 3:
672 label = "deltaR";
e7fa65db 673 nbins = 10;
55f64c2f 674 xmin = 0.;
e7fa65db 675 xmax = 5.0;
55f64c2f 676 break;
677
55f64c2f 678 case 4:
679 label = "deltaEta";
22d3f9c2 680 nbins = 24;
681 xmin = -1.2;
682 xmax = 1.2;
55f64c2f 683 break;
684
55f64c2f 685 case 5:
686 label = "deltaPhi";
9e4a0350 687 nbins = 72;
55f64c2f 688 xmin = -0.5*pi;
689 xmax = 1.5*pi;
f569a5d2 690 break;
55f64c2f 691
692 case 6:
693 label = "leading track";
694 nbins = 13;
695 xmin = 0;
696 xmax = 50;
697 break;
698
699 case 7:
700
701 label = "trigger track";
702 nbins =10;
703 xmin = 0;
704 xmax = 50;
705 break;
706
c1076e58 707 case 8:
708 label = "leading jet";
709 nbins = 3;
710 xmin = -0.5;
711 xmax = 2.5;
712 break;
55f64c2f 713
55f64c2f 714 }
55f64c2f 715}
716
55f64c2f 717//_________________________________________________
718// From CF event mixing code PhiCorrelations
719TObjArray* AliAnalysisTaskEmcalJetHMEC::CloneAndReduceTrackList(TObjArray* tracks)
720{
22d3f9c2 721 // clones a track list by using AliPicoTrack which uses much less memory (used for event mixing)
55f64c2f 722
723 TObjArray* tracksClone = new TObjArray;
724 tracksClone->SetOwner(kTRUE);
725
726 for (Int_t i=0; i<tracks->GetEntriesFast(); i++)
727 {
728 AliVParticle* particle = (AliVParticle*) tracks->At(i);
22d3f9c2 729 if(TMath::Abs(particle->Eta())>fTrkEta) continue;
f569a5d2 730 if(particle->Pt()<0.15) continue;
965c985f 731
22d3f9c2 732 Double_t trackpt=particle->Pt();
733
734 Int_t hadbin=-1;
735 if(trackpt<0.5) hadbin=0;
736 else if(trackpt<1) hadbin=1;
737 else if(trackpt<2) hadbin=2;
738 else if(trackpt<3) hadbin=3;
739 else if(trackpt<5) hadbin=4;
740 else if(trackpt<8) hadbin=5;
741 else if(trackpt<20) hadbin=6;
742
743 if(hadbin>-1) fHistTrackEtaPhi[hadbin]->Fill(particle->Eta(),particle->Phi());
744
5b38fc6d 745 tracksClone->Add(new AliPicoTrack(particle->Pt(), particle->Eta(), particle->Phi(), particle->Charge(), 0, 0, 0, 0));
55f64c2f 746 }
747
748 return tracksClone;
749}