Further update of the mu-h correlation task.
[u/mrichter/AliRoot.git] / PWGCF / Correlations / DPhi / MuonHadron / AliAnalysisTaskMuonHadronCorrelations.cxx
1 #include "AliLog.h"
2 #include "TH1D.h"
3 #include "TH2D.h"
4 #include "TTree.h"
5 #include "TAxis.h"
6 #include "AliAODEvent.h"
7 #include "AliAODTrack.h"
8 #include "TMath.h"
9 #include "TString.h"
10
11 #include "AliAnalysisTaskMuonHadronCorrelations.h"
12 #include "AliAnalysisManager.h"
13 #include "AliInputEventHandler.h"
14 #include "AliEventPoolManager.h"
15
16 ClassImp(AliAnalysisTaskMuonHadronCorrelations)
17
18 //====================================================================================================================================================
19
20 AliAnalysisTaskMuonHadronCorrelations::AliAnalysisTaskMuonHadronCorrelations() : 
21   AliAnalysisTaskSE(), 
22   fAOD(0x0),
23   fPoolMgr(0x0),
24   fTrackCB(0x0),
25   fTrackMA(0x0),
26   fFilterBitCentralBarrel(0),
27   fMaxEtaCentralBarrel(1.0),
28   fMaxChi2Muon(9999999999.), 
29   fMinRAbsMuon(0), 
30   fMaxRAbsMuon(9999999999.),
31   fTriggerMatchLevelMuon(0),
32   fNbinsCent(1), 
33   fNbinsPt(1),
34   fCentAxis(0x0), 
35   fPtAxis(0x0),
36   fHistV0Multiplicity(0x0), 
37   fHistITSMultiplicity(0x0),
38   fHistCentrality(0x0),
39   fHistEvStat(0x0),
40   fCentMethod(0),
41   fOutputList(0x0)
42 {
43
44   // Default constructor
45
46   for (Int_t iCent=0; iCent<fNMaxBinsCentrality; iCent++) {
47     for (Int_t iPtBinCB=0; iPtBinCB<fNMaxBinsPt; iPtBinCB++) {
48       for (Int_t iPtBinMA=0; iPtBinMA<fNMaxBinsPt; iPtBinMA++) {
49         fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA]    = NULL;
50         fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA] = NULL;
51       }
52     }
53     fHistNTracksCB_vs_NTracksMA[iCent]  = NULL;
54     fHistNTracksCB_vs_NTracksMAmixed[iCent]  = NULL;
55     fHistSingleMuonsPt[iCent]   = NULL;
56     fHistSingleMuonsPtmixed[iCent]   = NULL;
57     fHistSingleMuonsEtaVsPt[iCent]   = NULL;
58     fHistSingleMuonsEtaVsRAbs[iCent] = NULL;
59   }  
60   
61 }
62
63
64 //====================================================================================================================================================
65
66 AliAnalysisTaskMuonHadronCorrelations::AliAnalysisTaskMuonHadronCorrelations(const char *name) : 
67   AliAnalysisTaskSE(name), 
68   fAOD(0x0),
69   fPoolMgr(0x0),
70   fTrackCB(0x0),
71   fTrackMA(0x0),
72   fFilterBitCentralBarrel(0),
73   fMaxEtaCentralBarrel(1.0),
74   fMaxChi2Muon(9999999999.), 
75   fMinRAbsMuon(0), 
76   fMaxRAbsMuon(9999999999.),
77   fTriggerMatchLevelMuon(0),
78   fNbinsCent(1), 
79   fNbinsPt(1),
80   fCentAxis(0x0), 
81   fPtAxis(0x0),
82   fHistV0Multiplicity(0x0), 
83   fHistITSMultiplicity(0x0),
84   fHistCentrality(0x0),
85   fHistEvStat(0x0),
86   fCentMethod(0),
87   fOutputList(0x0)
88 {
89
90   // Constructor
91
92   for (Int_t iCent=0; iCent<fNMaxBinsCentrality; iCent++) {
93     for (Int_t iPtBinCB=0; iPtBinCB<fNMaxBinsPt; iPtBinCB++) {
94       for (Int_t iPtBinMA=0; iPtBinMA<fNMaxBinsPt; iPtBinMA++) {
95         fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA]    = NULL;
96         fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA] = NULL;
97       }
98     }
99     fHistNTracksCB_vs_NTracksMA[iCent]  = NULL;
100     fHistNTracksCB_vs_NTracksMAmixed[iCent]  = NULL;
101     fHistSingleMuonsPt[iCent]   = NULL;
102     fHistSingleMuonsPtmixed[iCent]   = NULL;
103     fHistSingleMuonsEtaVsPt[iCent]   = NULL;
104     fHistSingleMuonsEtaVsRAbs[iCent] = NULL;
105   }  
106   
107   // Define input and output slots here
108   DefineOutput(1, TList::Class());
109   
110 }
111
112 //====================================================================================================================================================
113
114 AliAnalysisTaskMuonHadronCorrelations::~AliAnalysisTaskMuonHadronCorrelations() {
115   
116   delete fCentAxis;
117   delete fPtAxis;
118
119   if (fOutputList  && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) 
120     delete fOutputList;
121 }
122
123 //====================================================================================================================================================
124
125 void AliAnalysisTaskMuonHadronCorrelations::UserCreateOutputObjects() {
126   
127   fOutputList = new TList();
128   fOutputList->SetOwner(kTRUE);
129
130   for (Int_t iCent=0; iCent<fNbinsCent; iCent++) {
131     for (Int_t iPtBinCB=0; iPtBinCB<fNbinsPt; iPtBinCB++) {
132       for (Int_t iPtBinMA=0; iPtBinMA<fNbinsPt; iPtBinMA++) {
133         
134         fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA]    = new TH1D(Form("fHistDeltaPhi_Cent%02d_PtBin%02d_%02d",iCent,iPtBinCB,iPtBinMA), 
135                                                                Form("%d-%d %%, %3.1f<p_{T}^{TPC}<%3.1f, %3.1f<p_{T}^{Muon}<%3.1f",
136                                                                     Int_t(fCentAxis->GetBinLowEdge(iCent+1)),
137                                                                     Int_t(fCentAxis->GetBinUpEdge(iCent+1)),
138                                                                     fPtAxis->GetBinLowEdge(iPtBinCB+1),
139                                                                     fPtAxis->GetBinUpEdge(iPtBinCB+1),
140                                                                     fPtAxis->GetBinLowEdge(iPtBinMA+1),
141                                                                     fPtAxis->GetBinUpEdge(iPtBinMA+1)),
142                                                                100, -0.5*TMath::RadToDeg()*TMath::Pi(), 1.5*TMath::RadToDeg()*TMath::Pi());
143         
144         fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA]    = new TH1D(Form("fHistDeltaPhiMix_Cent%02d_PtBin%02d_%02d",iCent,iPtBinCB,iPtBinMA), 
145                                                                   Form("%d-%d %%, %3.1f<p_{T}^{TPC}<%3.1f, %3.1f<p_{T}^{Muon}<%3.1f MIXED",
146                                                                        Int_t(fCentAxis->GetBinLowEdge(iCent+1)),
147                                                                        Int_t(fCentAxis->GetBinUpEdge(iCent+1)),
148                                                                        fPtAxis->GetBinLowEdge(iPtBinCB+1),
149                                                                        fPtAxis->GetBinUpEdge(iPtBinCB+1),
150                                                                        fPtAxis->GetBinLowEdge(iPtBinMA+1),
151                                                                        fPtAxis->GetBinUpEdge(iPtBinMA+1)),
152                                                                   100, -0.5*TMath::RadToDeg()*TMath::Pi(), 1.5*TMath::RadToDeg()*TMath::Pi());
153         
154         fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA]    -> SetXTitle("#Delta#varphi  [degrees]");
155         fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA] -> SetXTitle("#Delta#varphi  [degrees]");
156
157         fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA]    -> Sumw2();
158         fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA] -> Sumw2();
159         
160         fOutputList -> Add(fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA]);
161         fOutputList -> Add(fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA]);
162
163       }
164     }
165
166     fHistNTracksCB_vs_NTracksMA[iCent] = new TH2D(Form("fHistNTracksCB_vs_NTracksMA_Cent%02d",iCent),
167                                                   Form("%d-%d %%",Int_t(fCentAxis->GetBinLowEdge(iCent+1)),Int_t(fCentAxis->GetBinUpEdge(iCent+1))),
168                                                   100, 0, 500, 20, 0, 20);
169     fHistNTracksCB_vs_NTracksMA[iCent] -> SetXTitle("N_{tracks} Central Barrel");
170     fHistNTracksCB_vs_NTracksMA[iCent] -> SetYTitle("N_{tracks} Muon Arm");
171     fHistNTracksCB_vs_NTracksMA[iCent] -> Sumw2();
172
173     fHistNTracksCB_vs_NTracksMAmixed[iCent] = new TH2D(Form("fHistNTracksCB_vs_NTracksMAmixed_Cent%02d",iCent),
174                                                   Form("%d-%d %% MIXED",Int_t(fCentAxis->GetBinLowEdge(iCent+1)),Int_t(fCentAxis->GetBinUpEdge(iCent+1))),
175                                                   100, 0, 500, 20, 0, 20);
176     fHistNTracksCB_vs_NTracksMAmixed[iCent] -> SetXTitle("N_{tracks} Central Barrel");
177     fHistNTracksCB_vs_NTracksMAmixed[iCent] -> SetYTitle("N_{tracks} Muon Arm");
178     fHistNTracksCB_vs_NTracksMAmixed[iCent] -> Sumw2();
179
180     fOutputList -> Add(fHistNTracksCB_vs_NTracksMA[iCent]);
181     fOutputList -> Add(fHistNTracksCB_vs_NTracksMAmixed[iCent]);
182
183     fHistSingleMuonsPt[iCent] = new TH1D(Form("fHistSingleMuonPt_Cent%02d",iCent),
184                                          "p_{T} for single muons",
185                                          fNbinsPt, (Double_t*)fPtAxis->GetXbins()->GetArray());
186     fHistSingleMuonsPtmixed[iCent] = new TH1D(Form("fHistSingleMuonPtmixed_Cent%02d",iCent),
187                                               "p_{T} for single muons",
188                                               fNbinsPt, (Double_t*)fPtAxis->GetXbins()->GetArray());
189     fOutputList -> Add(fHistSingleMuonsPt[iCent]);
190     fOutputList -> Add(fHistSingleMuonsPtmixed[iCent]);
191
192     fHistSingleMuonsEtaVsPt[iCent]   = new TH2D(Form("fHistSingleMuonsEtaVsPt_Cent%02d",iCent),
193                                                "#eta vs p_{T} for single muons",
194                                                 100, -4.5, -2., 100, 0., 10.);
195     fOutputList->Add(fHistSingleMuonsEtaVsPt[iCent]);
196     fHistSingleMuonsEtaVsRAbs[iCent] = new TH2D(Form("fHistSingleMuonsEtaVsRAbs_Cent%02d",iCent),
197                                                 "#eta vs R_{Abs} for single muons",
198                                                 100, -4.5, -2., 100, 0, 100);
199     fOutputList->Add(fHistSingleMuonsEtaVsRAbs[iCent]);
200     
201   }
202   
203   fHistV0Multiplicity = new TH1D("fHistV0Multiplicity", "V0 Multiplicity", 500, 0, 1000);
204   fHistV0Multiplicity -> SetXTitle("Multiplicity");
205   fHistV0Multiplicity -> Sumw2();
206
207   fHistITSMultiplicity = new TH1D("fHistITSMultiplicity", "ITS Multiplicity", 500, 0, 500);
208   fHistITSMultiplicity -> SetXTitle("N_{Clusters}");
209   fHistITSMultiplicity -> Sumw2();
210
211   fHistCentrality = new TH1D("fHistCentrality", Form("%s Centrality",fCentMethod.Data()), 300, -100, 200);
212   fHistCentrality -> SetXTitle("Centrality  [%]");
213   fHistCentrality -> Sumw2();
214
215   fOutputList -> Add(fHistV0Multiplicity);
216   fOutputList -> Add(fHistITSMultiplicity);
217   fOutputList -> Add(fHistCentrality);
218
219   fHistEvStat = new TH1D("fHistEvStat","Event cuts statistics",20,-0.5,19.5);
220   fHistEvStat->SetXTitle("Cut index");
221   fOutputList->Add(fHistEvStat);
222
223   const Int_t kNZvtxBins  = 10;
224   // bins for further buffers are shifted by 100 cm
225   Double_t vertexBins[kNZvtxBins+1] = { -10,   -8,  -6,  -4,  -2,   0,   2,   4,   6,   8,  10 };
226   Int_t nZvtxBins  = kNZvtxBins;
227   Double_t* zvtxbin = vertexBins;
228
229   fPoolMgr = new AliEventPoolManager(1000, 20000, fNbinsCent, (Double_t*)fCentAxis->GetXbins()->GetArray(), nZvtxBins, zvtxbin);
230
231   PostData(1, fOutputList); 
232
233 }
234
235 //====================================================================================================================================================
236
237 void AliAnalysisTaskMuonHadronCorrelations::UserExec(Option_t *) {
238
239   fAOD = dynamic_cast<AliAODEvent *>(InputEvent());  
240   if (!fAOD) return;  
241
242   Int_t cutIndex = 0;
243   fHistEvStat->Fill(cutIndex++);
244   // Trigger selection
245   if (!(IsTriggerFired())) return;
246   fHistEvStat->Fill(cutIndex++);
247   
248   fHistV0Multiplicity  -> Fill(GetV0Multiplicity());
249   fHistITSMultiplicity -> Fill(GetITSMultiplicity());
250
251   Int_t centBin = GetCentBin();
252   if (centBin<0) return;
253   fHistEvStat->Fill(cutIndex++);
254   Double_t percentile = fAOD->GetCentrality()->GetCentralityPercentile(fCentMethod.Data());
255   fHistCentrality->Fill(percentile);
256
257   // Vertex selection
258   const AliAODVertex* trkVtx = fAOD->GetPrimaryVertex();
259   if (!trkVtx || trkVtx->GetNContributors()<=0) return;
260   TString vtxTtl = trkVtx->GetTitle();
261   if (!vtxTtl.Contains("VertexerTracks")) return;
262   fHistEvStat->Fill(cutIndex++);
263   Double_t zvtx = trkVtx->GetZ();
264   const AliAODVertex* spdVtx = fAOD->GetPrimaryVertexSPD();
265   if (spdVtx->GetNContributors()<=0) return;
266   TString vtxTyp = spdVtx->GetTitle();
267   Double_t cov[6]={0};
268   spdVtx->GetCovarianceMatrix(cov);
269   Double_t zRes = TMath::Sqrt(cov[5]);
270   if (vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) return;
271   if (TMath::Abs(spdVtx->GetZ() - trkVtx->GetZ())>0.5) return;
272   fHistEvStat->Fill(cutIndex++);
273
274   if (TMath::Abs(zvtx) > 10.) return;
275   fHistEvStat->Fill(cutIndex++);
276
277   TObjArray *tracksMuonArm = GetAcceptedTracksMuonArm(fAOD,centBin);
278   if (tracksMuonArm->GetEntriesFast() == 0) {
279     delete tracksMuonArm;
280     return;
281   }
282   fHistEvStat->Fill(cutIndex++);
283   TObjArray *tracksCentralBarrel = GetAcceptedTracksCentralBarrel(fAOD);
284   
285   fHistNTracksCB_vs_NTracksMA[centBin] -> Fill(tracksCentralBarrel->GetEntries(), tracksMuonArm->GetEntries());
286
287   AliDebug(1, Form("Single Event analysis : nTracksCB = %4d, nTracksMA = %4d", tracksCentralBarrel->GetEntries(), tracksMuonArm->GetEntries()));
288
289   // Same event  
290   for (Int_t iTrMA=0; iTrMA<tracksMuonArm->GetEntriesFast(); iTrMA++) {
291     fTrackMA = (AliAODTrack*) tracksMuonArm->At(iTrMA);
292     fHistSingleMuonsPt[centBin]->Fill(fTrackMA->Pt());
293     for (Int_t iTrCB=0; iTrCB<tracksCentralBarrel->GetEntriesFast(); iTrCB++) {
294       fTrackCB = (AliAODTrack*) tracksCentralBarrel -> At(iTrCB);
295       FillHistograms(centBin, kSingleEvent);
296     }
297   }
298
299   // Mixed event
300   {
301     AliEventPool* pool = fPoolMgr->GetEventPool(percentile, zvtx);
302     //pool->PrintInfo();
303     if (pool->IsReady() || pool->NTracksInPool() > 2000 || pool->GetCurrentNEvents() >= 5) 
304       for (Int_t jMix=0; jMix<pool->GetCurrentNEvents(); jMix++) {
305         TObjArray *mixedTracks = pool->GetEvent(jMix);
306         fHistNTracksCB_vs_NTracksMAmixed[centBin]->Fill(mixedTracks->GetEntriesFast(), tracksMuonArm->GetEntriesFast());
307         for (Int_t iTrMA=0; iTrMA<tracksMuonArm->GetEntriesFast(); iTrMA++) {
308           fTrackMA = (AliAODTrack*) tracksMuonArm->At(iTrMA);
309           fHistSingleMuonsPtmixed[centBin]->Fill(fTrackMA->Pt());
310           for (Int_t iTrCB=0; iTrCB<mixedTracks->GetEntriesFast(); iTrCB++) {
311             fTrackCB = (AliAODTrack*) mixedTracks -> At(iTrCB);
312             FillHistograms(centBin, kMixedEvent);
313           }
314         }
315       }
316     pool->UpdatePool(tracksCentralBarrel);
317   }
318
319   delete tracksMuonArm;
320
321   PostData(1, fOutputList); 
322
323 }
324
325 //====================================================================================================================================================
326
327 void AliAnalysisTaskMuonHadronCorrelations::FillHistograms(Int_t centrality, Int_t option) {
328
329   Int_t ptBinTrackCB = fPtAxis -> FindBin(fTrackCB->Pt());
330   Int_t ptBinTrackMA = fPtAxis -> FindBin(fTrackMA->Pt());
331
332   if (ptBinTrackCB<1 || ptBinTrackCB>fNbinsPt || ptBinTrackMA<1 || ptBinTrackMA>fNbinsPt) return;
333
334   Double_t deltaPhi = fTrackCB->Phi() - fTrackMA->Phi();
335   if (deltaPhi >  1.5*TMath::Pi()) deltaPhi -= TMath::TwoPi();
336   if (deltaPhi < -0.5*TMath::Pi()) deltaPhi += TMath::TwoPi();
337
338   if (option==kSingleEvent)     fHistDeltaPhi[centrality][ptBinTrackCB-1][ptBinTrackMA-1]    -> Fill(TMath::RadToDeg()*deltaPhi);
339   else if (option==kMixedEvent) fHistDeltaPhiMix[centrality][ptBinTrackCB-1][ptBinTrackMA-1] -> Fill(TMath::RadToDeg()*deltaPhi);
340
341 }
342
343 //====================================================================================================================================================
344
345 Bool_t AliAnalysisTaskMuonHadronCorrelations::IsTriggerFired() {
346   
347   Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kINT7); 
348
349   return isSelected;
350 }
351
352 //====================================================================================================================================================
353
354 Float_t AliAnalysisTaskMuonHadronCorrelations::GetV0Multiplicity() {
355   
356   Float_t multiplicity=0;
357   for (Int_t iChannel=0; iChannel<64; iChannel++) multiplicity += fAOD->GetVZEROData()->GetMultiplicity(iChannel);
358   return multiplicity;
359
360 }
361
362 //====================================================================================================================================================
363
364 TObjArray* AliAnalysisTaskMuonHadronCorrelations::GetAcceptedTracksCentralBarrel(AliAODEvent *aodEvent) {
365
366   // fills the array of central barrel tracks that pass the cuts
367
368   TObjArray *tracks = new TObjArray;
369   tracks->SetOwner(kTRUE);
370
371   Int_t nTracks = aodEvent->GetNTracks();
372
373   AliAODTrack *track = 0;
374
375   for (Int_t iTrack=0; iTrack<nTracks; iTrack++) {
376     track = aodEvent->GetTrack(iTrack);
377     if (track->TestFilterBit(fFilterBitCentralBarrel) && TMath::Abs(track->Eta())<fMaxEtaCentralBarrel) {
378       tracks->Add(new AliAODTrack(*track));
379     }
380   }
381
382   return tracks;
383
384 }
385
386 //====================================================================================================================================================
387
388 TObjArray* AliAnalysisTaskMuonHadronCorrelations::GetAcceptedTracksMuonArm(AliAODEvent *aodEvent, Int_t centBin) {
389
390   // fills the array of muon tracks that pass the cuts
391
392   TObjArray *tracks = new TObjArray;
393   tracks->SetOwner(kFALSE);
394
395   Int_t nTracks = aodEvent->GetNTracks();
396
397   AliAODTrack *track = 0;
398   
399   for (Int_t iTrack=0; iTrack<nTracks; iTrack++) {
400     track = aodEvent->GetTrack(iTrack);
401     if (track->IsMuonTrack() && track->GetMatchTrigger()>=fTriggerMatchLevelMuon) {
402       if (track->Chi2perNDF() < 5.) {
403         // histos
404         fHistSingleMuonsEtaVsPt[centBin]->Fill(track->Eta(),track->Pt());
405         fHistSingleMuonsEtaVsRAbs[centBin]->Fill(track->Eta(),track->GetRAtAbsorberEnd());
406         if (track->Eta() > -4. && track->Eta() < -2.5) {
407           Double_t rabs = track->GetRAtAbsorberEnd();
408           if (rabs > 17.6 && rabs < 89.5) {
409             tracks->Add(new AliAODTrack(*track));
410           }
411         }
412       }
413     }
414   }
415
416   return tracks;
417
418 }
419
420 //====================================================================================================================================================
421
422 void AliAnalysisTaskMuonHadronCorrelations::SetCentBinning(Int_t nBins, Double_t *limits) {
423
424   if (nBins>fNMaxBinsCentrality) {
425     AliInfo(Form("WARNING : only %d centrality bins (out of the %d proposed) will be considered",fNMaxBinsCentrality,nBins));
426     nBins = fNMaxBinsCentrality;
427   }
428   if (nBins<=0) {
429     AliInfo("WARNING : at least one centrality bin must be considered");
430     nBins = 1;
431   }
432   
433   fNbinsCent = nBins;
434   fCentAxis  = new TAxis(fNbinsCent, limits);
435
436 }
437
438 //====================================================================================================================================================
439
440 void AliAnalysisTaskMuonHadronCorrelations::SetPtBinning(Int_t nBins, Double_t *limits) {
441
442   if (nBins>fNMaxBinsPt) {
443     AliInfo(Form("WARNING : only %d pt bins (out of the %d proposed) will be considered",fNMaxBinsPt,nBins));
444     nBins = fNMaxBinsPt;
445   }
446   if (nBins<=0) {
447     AliInfo("WARNING : at least one pt bin must be considered");
448     nBins = 1;
449   }
450   
451   fNbinsPt = nBins;
452   fPtAxis  = new TAxis(fNbinsPt, limits);
453
454 }
455
456 //====================================================================================================================================================
457
458 Int_t AliAnalysisTaskMuonHadronCorrelations::GetCentBin() {
459
460   Double_t percentile = fAOD->GetCentrality()->GetCentralityPercentile(fCentMethod.Data());
461
462   Int_t bin = fCentAxis->FindBin(percentile) - 1;
463   if (bin >= fNbinsCent) bin = -1;
464   return bin;
465   
466 }
467
468 //====================================================================================================================================================
469
470 Double_t AliAnalysisTaskMuonHadronCorrelations::GetITSMultiplicity() {
471
472   Double_t multiplicity = fAOD->GetHeader()->GetNumberOfITSClusters(1);
473
474   return multiplicity;
475
476 }
477
478 //====================================================================================================================================================
479
480 void AliAnalysisTaskMuonHadronCorrelations::Terminate(Option_t *) {
481
482   fOutputList = dynamic_cast<TList*> (GetOutputData(1));
483   if (!fOutputList) {
484     printf("ERROR: Output list not available\n");
485     return;
486   }  
487
488 }
489
490 //====================================================================================================================================================
491