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