]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetHMEC.cxx
Remove hard-coded number of centrality bins
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskEmcalJetHMEC.cxx
... / ...
CommitLineData
1// $Id$
2
3//////////
4//Measure Jet-hadron correlations
5//Does event Mixing using AliEventPoolManager
6/////////
7
8#include "AliAnalysisTaskEmcalJetHMEC.h"
9
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"
34
35#include "TVector3.h"
36#include "AliPicoTrack.h"
37#include "AliEventPoolManager.h"
38
39ClassImp(AliAnalysisTaskEmcalJetHMEC)
40
41//________________________________________________________________________
42AliAnalysisTaskEmcalJetHMEC::AliAnalysisTaskEmcalJetHMEC() :
43 AliAnalysisTaskEmcalJet("HMEC",kFALSE),
44 fTracksName(""),
45 fJetsName(""),
46 fPhimin(-10),
47 fPhimax(10),
48 fEtamin(-0.9),
49 fEtamax(0.9),
50 fAreacut(0.0),
51 fTrkBias(5),
52 fClusBias(5),
53 fTrkEta(0.9),
54 fDoEventMixing(0),
55 fMixingTracks(50000),
56 fESD(0),
57 fAOD(0),
58 fPoolMgr(0x0),
59 fHistTrackPt(0),
60 fHistCentrality(0),
61 fHistJetEtaPhi(0),
62 fHistJetHEtaPhi(0),
63 fhnMixedEvents(0x0),
64 fhnJH(0x0)
65{
66 // Default Constructor
67
68 for(Int_t ipta=0; ipta<7; ipta++){
69 fHistTrackEtaPhi[ipta]=0;
70 }
71
72 for(Int_t icent = 0; icent<6; ++icent){
73 fHistJetPt[icent]=0;
74 fHistJetPtBias[icent]=0;
75 fHistLeadJetPt[icent]=0;
76 fHistLeadJetPtBias[icent]=0;
77 fHistJetPtTT[icent]=0;
78 for(Int_t iptjet = 0; iptjet<5; ++iptjet){
79 for(Int_t ieta = 0; ieta<3; ++ieta){
80 fHistJetH[icent][iptjet][ieta]=0;
81 fHistJetHBias[icent][iptjet][ieta]=0;
82 fHistJetHTT[icent][iptjet][ieta]=0;
83 }
84 }
85 }
86
87}
88//________________________________________________________________________
89AliAnalysisTaskEmcalJetHMEC::AliAnalysisTaskEmcalJetHMEC(const char *name) :
90 AliAnalysisTaskEmcalJet(name,kTRUE),
91 fTracksName(""),
92 fJetsName(""),
93 fPhimin(-10),
94 fPhimax(10),
95 fEtamin(-0.9),
96 fEtamax(0.9),
97 fAreacut(0.0),
98 fTrkBias(5),
99 fClusBias(5),
100 fTrkEta(0.9),
101 fDoEventMixing(0),
102 fMixingTracks(50000),
103 fESD(0),
104 fAOD(0),
105 fPoolMgr(0x0),
106 fHistTrackPt(0),
107 fHistCentrality(0),
108 fHistJetEtaPhi(0),
109 fHistJetHEtaPhi(0),
110 fhnMixedEvents(0x0),
111 fhnJH(0x0)
112{
113 // Constructor
114 for(Int_t ipta=0; ipta<7; ipta++){
115 fHistTrackEtaPhi[ipta]=0;
116 }
117 for(Int_t icent = 0; icent<6; ++icent){
118 fHistJetPt[icent]=0;
119 fHistJetPtBias[icent]=0;
120 fHistLeadJetPt[icent]=0;
121 fHistLeadJetPtBias[icent]=0;
122 fHistJetPtTT[icent]=0;
123 for(Int_t iptjet = 0; iptjet<5; ++iptjet){
124 for(Int_t ieta = 0; ieta<3; ++ieta){
125 fHistJetH[icent][iptjet][ieta]=0;
126 fHistJetHBias[icent][iptjet][ieta]=0;
127 fHistJetHTT[icent][iptjet][ieta]=0;
128 }
129 }
130 }
131
132}
133
134//________________________________________________________________________
135void AliAnalysisTaskEmcalJetHMEC::UserCreateOutputObjects() {
136 // Called once
137 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
138 OpenFile(1);
139
140 // Create histograms
141 fHistTrackPt = new TH1F("fHistTrackPt", "P_{T} distribution", 1000, 0.0, 100.0);
142 fHistCentrality = new TH1F("fHistCentrality","centrality",100,0,100);
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);
145
146 TString name;
147
148 for(Int_t ipta=0; ipta<7; ++ipta){
149 name = Form("fHistTrackEtaPhi_%i", ipta);
150 fHistTrackEtaPhi[ipta] = new TH2F(name,name,400,-1,1,720,0.0,2.0*TMath::Pi());
151 fOutput->Add(fHistTrackEtaPhi[ipta]);
152 }
153
154 for(Int_t icent = 0; icent<6; ++icent){
155 name = Form("fHistJetPt_%i",icent);
156 fHistJetPt[icent] = new TH1F(name,name,200,0,200);
157 fOutput->Add(fHistJetPt[icent]);
158
159 name = Form("fHistJetPtBias_%i",icent);
160 fHistJetPtBias[icent] = new TH1F(name,name,200,0,200);
161 fOutput->Add(fHistJetPtBias[icent]);
162
163 name = Form("fHistLeadJetPt_%i",icent);
164 fHistLeadJetPt[icent] = new TH1F(name,name,200,0,200);
165 fOutput->Add(fHistLeadJetPt[icent]);
166
167 name = Form("fHistLeadJetPtBias_%i",icent);
168 fHistLeadJetPtBias[icent] = new TH1F(name,name,200,0,200);
169 fOutput->Add(fHistLeadJetPtBias[icent]);
170
171 name = Form("fHistJetPtTT_%i",icent);
172 fHistJetPtTT[icent] = new TH1F(name,name,200,0,200);
173 fOutput->Add(fHistJetPtTT[icent]);
174
175 for(Int_t iptjet = 0; iptjet<5; ++iptjet){
176 for(Int_t ieta = 0; ieta<3; ++ieta){
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]);
180
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]);
184
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]);
188
189 }
190 }
191 }
192
193 UInt_t cifras = 0; // bit coded, see GetDimParams() below
194 cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7 | 1<<8;
195 fhnJH = NewTHnSparseF("fhnJH", cifras);
196 fhnJH->Sumw2();
197 fOutput->Add(fhnJH);
198
199 if(fDoEventMixing){
200 cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7 | 1<<8;
201 fhnMixedEvents = NewTHnSparseF("fhnMixedEvents", cifras);
202 fhnMixedEvents->Sumw2();
203 fOutput->Add(fhnMixedEvents);
204 }
205
206 fOutput->Add(fHistTrackPt);
207 fOutput->Add(fHistCentrality);
208 fOutput->Add(fHistJetEtaPhi);
209 fOutput->Add(fHistJetHEtaPhi);
210
211 PostData(1, fOutput);
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
217 Int_t nZvtxBins = 5+1+5;
218 // bins for second buffer are shifted by 100 cm
219 Double_t vertexBins[] = { -10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10, };
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 }
227
228 fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBins, centralityBins, nZvtxBins, zvtxbin);
229
230}
231
232//________________________________________________________________________
233
234Double_t AliAnalysisTaskEmcalJetHMEC:: RelativePhi(Double_t mphi,Double_t vphi) {
235
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
247//________________________________________________________________________
248Int_t AliAnalysisTaskEmcalJetHMEC::GetCentBin(Double_t cent) const {
249 // Get centrality bin.
250
251 Int_t centbin = -1;
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;
258 return centbin;
259}
260
261//________________________________________________________________________
262Int_t AliAnalysisTaskEmcalJetHMEC::GetEtaBin(Double_t eta) const {
263 // Get eta bin for histos.
264
265 Int_t etabin = -1;
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;
269 return etabin;
270}
271//________________________________________________________________________
272Int_t AliAnalysisTaskEmcalJetHMEC::GetpTjetBin(Double_t pt) const
273{
274 // Get jet pt bin for histos.
275
276 Int_t ptbin = -1;
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;
282
283 return ptbin;
284}
285
286//________________________________________________________________________
287void AliAnalysisTaskEmcalJetHMEC::ExecOnce() {
288 AliAnalysisTaskEmcalJet::ExecOnce();
289
290}
291
292//________________________________________________________________________
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 }
303
304 // what kind of event do we have: AOD or ESD?
305 Bool_t esdMode = kTRUE;
306 if (dynamic_cast<AliAODEvent*>(InputEvent())) esdMode = kFALSE;
307
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 }
316
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 }
325
326 TList *list = InputEvent()->GetList();
327 if(!list) {
328 AliError(Form("ERROR: list not attached\n"));
329 return kTRUE;
330 }
331
332 // get centrality
333 if (fCent<0) {
334 AliError(Form("Centrality negative: %f", fCent));
335 return kTRUE;
336 }
337
338 Int_t centbin = GetCentBin(fCent);
339 if(centbin<0) return kTRUE;
340
341 Double_t fvertex[3]={0,0,0};
342 InputEvent()->GetPrimaryVertex()->GetXYZ(fvertex);
343 Double_t zVtx=fvertex[2];
344
345 if(fabs(zVtx)>10.0) return kTRUE;
346
347 fHistCentrality->Fill(fCent);
348
349 TClonesArray *jets = 0;
350 TClonesArray *tracks = 0;
351
352 tracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracks));
353 if (!tracks) {
354 AliError(Form("Pointer to tracks %s == 0", fTracks->GetName() ));
355 return kTRUE;
356 }
357 const Int_t Ntracks=tracks->GetEntries();
358
359 jets= dynamic_cast<TClonesArray*>(list->FindObject(fJets));
360 if (!jets) {
361 AliError(Form("Pointer to tracks %s == 0", fJets->GetName() ));
362 return kTRUE;
363 }
364 const Int_t Njets = jets->GetEntries();
365
366 //Leticia's loop to find hardest track
367 Int_t iTT=-1;
368 Double_t ptmax=-10;
369
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 }
376
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;
383 }
384 }
385
386 Int_t ijethi=-1;
387 Double_t highestjetpt=0.0;
388 Int_t passedTTcut=0;
389
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;
394
395 Double_t jetPt = jet->Pt();
396
397 if(highestjetpt<jetPt){
398 ijethi=ijets;
399 highestjetpt=jetPt;
400 }
401 }
402
403 for (Int_t ijet = 0; ijet < Njets; ijet++){
404 AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet));
405 if (!jet) continue;
406 if(!AcceptthisJet(jet)) continue;
407
408 Double_t jetphi = jet->Phi();
409 Double_t jetPt = jet->Pt();
410 Double_t jeteta=jet->Eta();
411
412 Double_t leadjet=0;
413 if (ijet==ijethi) leadjet=1;
414
415 fHistJetPt[centbin]->Fill(jet->Pt());
416 fHistLeadJetPt[centbin]->Fill(jet->Pt());
417
418 if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){
419 fHistJetPtBias[centbin]->Fill(jet->Pt());
420 fHistLeadJetPtBias[centbin]->Fill(jet->Pt());
421 }
422
423 fHistJetEtaPhi->Fill(jet->Eta(),jetphi);
424
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 }
430
431 if(passedTTcut)
432 fHistJetPtTT[centbin]->Fill(jet->Pt());
433
434 Int_t iptjet=-1;
435 iptjet=GetpTjetBin(jetPt);
436 if(iptjet<0) continue;
437
438 //if (highestjetpt>15) {
439 if (highestjetpt>8) {
440
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 }
447
448 if(TMath::Abs(track->Eta())>fTrkEta) continue;
449
450 fHistTrackPt->Fill(track->Pt());
451
452 if (track->Pt()<0.15) continue;
453
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 }
466
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();
473
474 fHistJetH[centbin][iptjet][ieta]->Fill(dphijh,track->Pt());
475 fHistJetHEtaPhi->Fill(deta,dphijh);
476
477 Double_t dR=sqrt(deta*deta+dphijh*dphijh);
478
479 if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){
480 fHistJetHBias[centbin][iptjet][ieta]->Fill(dphijh,trackpt);
481
482 Double_t triggerEntries[8] = {fCent,jetPt,trackpt,dR,deta,dphijh,0.0,leadjet};
483 fhnJH->Fill(triggerEntries);
484 }
485
486 if(passedTTcut)
487 fHistJetHTT[centbin][iptjet][ieta]->Fill(dphijh,trackpt);
488
489 } //track loop
490 }//jet pt cut
491 }//jet loop
492
493 //Prepare to do event mixing
494
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
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
519 AliEventPool* pool = fPoolMgr->GetEventPool(fCent, zVtx);
520
521 if (!pool){
522 AliFatal(Form("No pool found for centrality = %f, zVtx = %f", fCent, zVtx));
523 return kTRUE;
524 }
525
526 //check for a trigger jet
527 if (pool->IsReady() || pool->NTracksInPool() > fMixingTracks / 10 || pool->GetCurrentNEvents() >= 5) {
528
529 for (Int_t ijet = 0; ijet < Njets; ijet++){
530 Double_t leadjet=0;
531 if (ijet==ijethi) leadjet=1;
532
533 AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet));
534 if(!AcceptthisJet(jet)) continue;
535
536 Double_t jetPt = jet->Pt();
537 Double_t jetphi = jet->Phi();
538 Double_t jeteta=jet->Eta();
539
540 Int_t nMix = pool->GetCurrentNEvents();
541
542 //Fill for biased jet triggers only
543 if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){
544
545 // Fill mixed-event histos here
546 for (Int_t jMix=0; jMix<nMix; jMix++) {
547 TObjArray* bgTracks = pool->GetEvent(jMix);
548 const Int_t Nbgtrks = bgTracks->GetEntries();
549
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 }
563 }
564 }
565 }
566 }
567
568 //update pool if jet in event or not
569 pool->UpdatePool(tracksClone);
570 }
571
572 return kTRUE;
573}
574
575//________________________________________________________________________
576void AliAnalysisTaskEmcalJetHMEC::Terminate(Option_t *)
577{
578 //just terminate
579}
580
581//________________________________________________________________________
582Int_t AliAnalysisTaskEmcalJetHMEC::AcceptthisJet(AliEmcalJet *jet)
583{
584 //applies all jet cuts except pt
585 float jetphi = jet->Phi();
586 if (jetphi>TMath::Pi())
587 jetphi = jetphi-2*TMath::Pi();
588
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)
597 return 0;
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;
604}
605
606//________________________________________________________________________
607
608THnSparse* AliAnalysisTaskEmcalJetHMEC::NewTHnSparseF(const char* name, UInt_t entries){
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 (%)";
652 nbins = 10;
653 xmin = 0.;
654 xmax = 100.;
655 break;
656
657 case 1:
658 label = "corrected jet pt";
659 nbins = 20;
660 xmin = 0.;
661 xmax = 200.;
662 break;
663
664 case 2:
665 label = "track pT";
666 nbins = 100;
667 xmin = 0.;
668 xmax = 10;
669 break;
670
671 case 3:
672 label = "deltaR";
673 nbins = 10;
674 xmin = 0.;
675 xmax = 5.0;
676 break;
677
678 case 4:
679 label = "deltaEta";
680 nbins = 24;
681 xmin = -1.2;
682 xmax = 1.2;
683 break;
684
685 case 5:
686 label = "deltaPhi";
687 nbins = 72;
688 xmin = -0.5*pi;
689 xmax = 1.5*pi;
690 break;
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
707 case 8:
708 label = "leading jet";
709 nbins = 3;
710 xmin = -0.5;
711 xmax = 2.5;
712 break;
713
714 }
715}
716
717//_________________________________________________
718// From CF event mixing code PhiCorrelations
719TObjArray* AliAnalysisTaskEmcalJetHMEC::CloneAndReduceTrackList(TObjArray* tracks)
720{
721 // clones a track list by using AliPicoTrack which uses much less memory (used for event mixing)
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);
729 if(TMath::Abs(particle->Eta())>fTrkEta) continue;
730 if(particle->Pt()<0.15) continue;
731
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
745 tracksClone->Add(new AliPicoTrack(particle->Pt(), particle->Eta(), particle->Phi(), particle->Charge(), 0, 0, 0, 0));
746 }
747
748 return tracksClone;
749}