From 7308e5d4b0362bf74501e857da4aaa3ec3e82d4f Mon Sep 17 00:00:00 2001 From: cvetan Date: Thu, 7 Mar 2013 12:37:28 +0000 Subject: [PATCH 1/1] Adding muon-muon correlations task (Antonio) --- PWGCF/CMakelibPWGCFCorrelationsDPhi.pkg | 1 + .../AddAnalysisTaskDiMuonCorrelations.C | 41 ++ .../AliAnalysisTaskDiMuonCorrelations.cxx | 527 ++++++++++++++++++ .../AliAnalysisTaskDiMuonCorrelations.h | 100 ++++ PWGCF/PWGCFCorrelationsDPhiLinkDef.h | 1 + 5 files changed, 670 insertions(+) create mode 100644 PWGCF/Correlations/DPhi/MuonHadron/AddAnalysisTaskDiMuonCorrelations.C create mode 100644 PWGCF/Correlations/DPhi/MuonHadron/AliAnalysisTaskDiMuonCorrelations.cxx create mode 100644 PWGCF/Correlations/DPhi/MuonHadron/AliAnalysisTaskDiMuonCorrelations.h diff --git a/PWGCF/CMakelibPWGCFCorrelationsDPhi.pkg b/PWGCF/CMakelibPWGCFCorrelationsDPhi.pkg index 7defb21e191..47c0101a665 100644 --- a/PWGCF/CMakelibPWGCFCorrelationsDPhi.pkg +++ b/PWGCF/CMakelibPWGCFCorrelationsDPhi.pkg @@ -41,6 +41,7 @@ set ( SRCS Correlations/DPhi/AliLeadingV0Correlation.cxx Correlations/DPhi/AliAnalysisTaskLongRangeCorrelations.cxx Correlations/DPhi/MuonHadron/AliAnalysisTaskMuonHadronCorrelations.cxx + Correlations/DPhi/MuonHadron/AliAnalysisTaskDiMuonCorrelations.cxx ) string ( REPLACE ".cxx" ".h" HDRS "${SRCS}" ) diff --git a/PWGCF/Correlations/DPhi/MuonHadron/AddAnalysisTaskDiMuonCorrelations.C b/PWGCF/Correlations/DPhi/MuonHadron/AddAnalysisTaskDiMuonCorrelations.C new file mode 100644 index 00000000000..bf660082a15 --- /dev/null +++ b/PWGCF/Correlations/DPhi/MuonHadron/AddAnalysisTaskDiMuonCorrelations.C @@ -0,0 +1,41 @@ +AliAnalysisTaskDiMuonCorrelations *AddAnalysisTaskDiMuonCorrelations(const char *centMethod = "V0M") { + + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + if (!mgr) { + printf("Error in adding AnalysisTaskDiMuonCorrelations: no Analysis Manager found!\n"); + return NULL; + } + + AliAnalysisTaskDiMuonCorrelations *task = new AliAnalysisTaskDiMuonCorrelations(Form("AliAnalysisTaskDiMuonCorrelations_%s",centMethod)); + + // Set analysis cuts + task->SetEtaRangeMuon(-4.0, -2.5); + task->SetTriggerMatchLevelMuon(1); + + const Int_t nBinCent = 4; + Double_t centLimits[nBinCent+1] = {0., 20., 40, 60., 100.}; + task->SetCentBinning(nBinCent, centLimits); + + task->SetCentMethod(centMethod); + + const Int_t nBinPt = 3; + Double_t ptLimits[nBinPt+1] = {0., 1., 2., 4.}; + task->SetPtBinning(nBinPt, ptLimits); + + const Int_t nBinEta = 3; + Double_t etaLimits[nBinEta+1] = {-4., -3.6, -3.2, -2.5}; + task->SetEtaBinning(nBinEta, etaLimits); + + mgr->AddTask(task); + + // create output container + AliAnalysisDataContainer *output = mgr->CreateContainer(Form("DiMuonCorrHistos_%s",centMethod), TList::Class(), AliAnalysisManager::kOutputContainer, + Form("%s:DiMuonCorrelations_%s", AliAnalysisManager::GetCommonFileName(), centMethod)); + + // finaly connect input and output + mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer()); + mgr->ConnectOutput(task, 1, output); + + return task; +} + diff --git a/PWGCF/Correlations/DPhi/MuonHadron/AliAnalysisTaskDiMuonCorrelations.cxx b/PWGCF/Correlations/DPhi/MuonHadron/AliAnalysisTaskDiMuonCorrelations.cxx new file mode 100644 index 00000000000..2c82ba5b1b1 --- /dev/null +++ b/PWGCF/Correlations/DPhi/MuonHadron/AliAnalysisTaskDiMuonCorrelations.cxx @@ -0,0 +1,527 @@ +#include "AliLog.h" +#include "TH1D.h" +#include "TH2D.h" +#include "TTree.h" +#include "TAxis.h" +#include "AliAODEvent.h" +#include "AliAODTrack.h" +#include "TMath.h" +#include "TString.h" + +#include "AliAnalysisTaskDiMuonCorrelations.h" +#include "AliAnalysisManager.h" +#include "AliInputEventHandler.h" +#include "AliEventPoolManager.h" + +ClassImp(AliAnalysisTaskDiMuonCorrelations) + +//==================================================================================================================================================== + +AliAnalysisTaskDiMuonCorrelations::AliAnalysisTaskDiMuonCorrelations() : + AliAnalysisTaskSE(), + fAOD(0x0), + fPoolMgr(0x0), + fMaxChi2Muon(5.), + fMinEtaMuon(-4.0), + fMaxEtaMuon(-2.5), + fTriggerMatchLevelMuon(0), + fNbinsCent(1), + fNbinsPt(1), + fCentAxis(0x0), + fPtAxis(0x0), + fEtaAxis(0x0), + fHistV0Multiplicity(0x0), + fHistITSMultiplicity(0x0), + fHistCentrality(0x0), + fHistEvStat(0x0), + fCentMethod(0), + fOutputList(0x0) +{ + + // Default constructor + + fMuonTrack[0] = NULL; + fMuonTrack[1] = NULL; + + for (Int_t iCent=0; iCentIsProofMode()) delete fOutputList; + +} + +//==================================================================================================================================================== + +void AliAnalysisTaskDiMuonCorrelations::UserCreateOutputObjects() { + + fOutputList = new TList(); + fOutputList->SetOwner(kTRUE); + + for (Int_t iCent=0; iCentGetBinLowEdge(iCent+1)), + Int_t(fCentAxis->GetBinUpEdge(iCent+1)), + fPtAxis->GetBinLowEdge(iPtBinMuon1+1), + fPtAxis->GetBinUpEdge(iPtBinMuon1+1), + fPtAxis->GetBinLowEdge(iPtBinMuon2+1), + fPtAxis->GetBinUpEdge(iPtBinMuon2+1)), + 100, -0.5*TMath::RadToDeg()*TMath::Pi(), 1.5*TMath::RadToDeg()*TMath::Pi()); + + fHistDeltaPhiMix[iCent][iPtBinMuon1][iPtBinMuon2] = new TH1D(Form("fHistDeltaPhiMix_Cent%02d_PtBin%02d_%02d",iCent,iPtBinMuon1,iPtBinMuon2), + Form("%d-%d %%, %3.1fGetBinLowEdge(iCent+1)), + Int_t(fCentAxis->GetBinUpEdge(iCent+1)), + fPtAxis->GetBinLowEdge(iPtBinMuon1+1), + fPtAxis->GetBinUpEdge(iPtBinMuon1+1), + fPtAxis->GetBinLowEdge(iPtBinMuon2+1), + fPtAxis->GetBinUpEdge(iPtBinMuon2+1)), + 100, -0.5*TMath::RadToDeg()*TMath::Pi(), 1.5*TMath::RadToDeg()*TMath::Pi()); + + fHistDeltaPhi[iCent][iPtBinMuon1][iPtBinMuon2] -> SetXTitle("#Delta#varphi [degrees]"); + fHistDeltaPhiMix[iCent][iPtBinMuon1][iPtBinMuon2] -> SetXTitle("#Delta#varphi [degrees]"); + + fHistDeltaPhi[iCent][iPtBinMuon1][iPtBinMuon2] -> Sumw2(); + fHistDeltaPhiMix[iCent][iPtBinMuon1][iPtBinMuon2] -> Sumw2(); + + fOutputList -> Add(fHistDeltaPhi[iCent][iPtBinMuon1][iPtBinMuon2]); + fOutputList -> Add(fHistDeltaPhiMix[iCent][iPtBinMuon1][iPtBinMuon2]); + + fHistEtaDeltaPhi[iCent][iPtBinMuon1][iPtBinMuon2] = new TH2D(Form("fHistEtaDeltaPhi_Cent%02d_PtBin%02d_%02d",iCent,iPtBinMuon1,iPtBinMuon2), + Form("%d-%d %%, %3.1fGetBinLowEdge(iCent+1)), + Int_t(fCentAxis->GetBinUpEdge(iCent+1)), + fPtAxis->GetBinLowEdge(iPtBinMuon1+1), + fPtAxis->GetBinUpEdge(iPtBinMuon1+1), + fPtAxis->GetBinLowEdge(iPtBinMuon2+1), + fPtAxis->GetBinUpEdge(iPtBinMuon2+1)), + 100, -0.5*TMath::RadToDeg()*TMath::Pi(), 1.5*TMath::RadToDeg()*TMath::Pi(), + fEtaAxis->GetNbins(),(Double_t*)fEtaAxis->GetXbins()->GetArray()); + + fHistEtaDeltaPhiMix[iCent][iPtBinMuon1][iPtBinMuon2] = new TH2D(Form("fHistEtaDeltaPhiMix_Cent%02d_PtBin%02d_%02d",iCent,iPtBinMuon1,iPtBinMuon2), + Form("%d-%d %%, %3.1fGetBinLowEdge(iCent+1)), + Int_t(fCentAxis->GetBinUpEdge(iCent+1)), + fPtAxis->GetBinLowEdge(iPtBinMuon1+1), + fPtAxis->GetBinUpEdge(iPtBinMuon1+1), + fPtAxis->GetBinLowEdge(iPtBinMuon2+1), + fPtAxis->GetBinUpEdge(iPtBinMuon2+1)), + 100, -0.5*TMath::RadToDeg()*TMath::Pi(), 1.5*TMath::RadToDeg()*TMath::Pi(), + fEtaAxis->GetNbins(),(Double_t*)fEtaAxis->GetXbins()->GetArray()); + + fHistEtaDeltaPhi[iCent][iPtBinMuon1][iPtBinMuon2] -> SetXTitle("#Delta#varphi [degrees]"); + fHistEtaDeltaPhiMix[iCent][iPtBinMuon1][iPtBinMuon2] -> SetXTitle("#Delta#varphi [degrees]"); + fHistEtaDeltaPhi[iCent][iPtBinMuon1][iPtBinMuon2] -> SetYTitle("#eta"); + fHistEtaDeltaPhiMix[iCent][iPtBinMuon1][iPtBinMuon2] -> SetYTitle("#eta"); + + fHistEtaDeltaPhi[iCent][iPtBinMuon1][iPtBinMuon2] -> Sumw2(); + fHistEtaDeltaPhiMix[iCent][iPtBinMuon1][iPtBinMuon2] -> Sumw2(); + + fOutputList -> Add(fHistEtaDeltaPhi[iCent][iPtBinMuon1][iPtBinMuon2]); + fOutputList -> Add(fHistEtaDeltaPhiMix[iCent][iPtBinMuon1][iPtBinMuon2]); + + } + } + + fHistNMuons_vs_NMuons[iCent] = new TH2D(Form("fHistNMuons_vs_NMuons_Cent%02d",iCent), + Form("%d-%d %%",Int_t(fCentAxis->GetBinLowEdge(iCent+1)),Int_t(fCentAxis->GetBinUpEdge(iCent+1))), + 100, 0, 500, 20, 0, 20); + fHistNMuons_vs_NMuons[iCent] -> SetXTitle("N_{tracks} Muon Arm"); + fHistNMuons_vs_NMuons[iCent] -> SetYTitle("N_{tracks} Muon Arm"); + fHistNMuons_vs_NMuons[iCent] -> Sumw2(); + + fHistNMuons_vs_NMuons_Mixed[iCent] = new TH2D(Form("fHistNMuons_vs_NMuons_Mixed_Cent%02d",iCent), + Form("%d-%d %% MIXED",Int_t(fCentAxis->GetBinLowEdge(iCent+1)),Int_t(fCentAxis->GetBinUpEdge(iCent+1))), + 100, 0, 500, 20, 0, 20); + fHistNMuons_vs_NMuons_Mixed[iCent] -> SetXTitle("N_{tracks} Muon Arm"); + fHistNMuons_vs_NMuons_Mixed[iCent] -> SetYTitle("N_{tracks} Muon Arm"); + fHistNMuons_vs_NMuons_Mixed[iCent] -> Sumw2(); + + fHistTracksEtavsEta[iCent] = new TH2D(Form("fHistTracksEtavsEta_%02d",iCent), "#eta 1st muon vs #eta 2nd muon", 100,-4.5,-2.,100,-1.5,1.5); + fHistTracksEtavsEta_Mixed[iCent] = new TH2D(Form("fHistTracksEtavsEta_Mixed_%02d",iCent), "#eta 1st muon vs #eta 2nd muon", 100,-4.5,-2.,100,-1.5,1.5); + + fOutputList -> Add(fHistNMuons_vs_NMuons[iCent]); + fOutputList -> Add(fHistNMuons_vs_NMuons_Mixed[iCent]); + fOutputList -> Add(fHistTracksEtavsEta[iCent]); + fOutputList -> Add(fHistTracksEtavsEta_Mixed[iCent]); + + fHistSingleMuonsPt[iCent] = new TH1D(Form("fHistSingleMuonPt_Cent%02d",iCent), "p_{T} for single muons", fNbinsPt, (Double_t*)fPtAxis->GetXbins()->GetArray()); + fHistSingleMuonsPt_Mixed[iCent] = new TH1D(Form("fHistSingleMuonPtmixed_Cent%02d",iCent), "p_{T} for single muons", fNbinsPt, (Double_t*)fPtAxis->GetXbins()->GetArray()); + fOutputList -> Add(fHistSingleMuonsPt[iCent]); + fOutputList -> Add(fHistSingleMuonsPt_Mixed[iCent]); + + fHistSingleMuonsEtaPt[iCent] = new TH2D(Form("fHistSingleMuonEtaPt_Cent%02d",iCent), "#eta vs p_{T} for single muons", + fNbinsPt, (Double_t*)fPtAxis->GetXbins()->GetArray(), + fEtaAxis->GetNbins(),(Double_t*)fEtaAxis->GetXbins()->GetArray()); + fHistSingleMuonsEtaPt_Mixed[iCent] = new TH2D(Form("fHistSingleMuonEtaPtmixed_Cent%02d",iCent), "#eta vs p_{T} for single muons", + fNbinsPt, (Double_t*)fPtAxis->GetXbins()->GetArray(), + fEtaAxis->GetNbins(),(Double_t*)fEtaAxis->GetXbins()->GetArray()); + fOutputList -> Add(fHistSingleMuonsEtaPt[iCent]); + fOutputList -> Add(fHistSingleMuonsEtaPt_Mixed[iCent]); + + } + + fHistV0Multiplicity = new TH1D("fHistV0Multiplicity", "V0 Multiplicity", 500, 0, 1000); + fHistV0Multiplicity -> SetXTitle("Multiplicity"); + fHistV0Multiplicity -> Sumw2(); + + fHistITSMultiplicity = new TH1D("fHistITSMultiplicity", "ITS Multiplicity", 500, 0, 500); + fHistITSMultiplicity -> SetXTitle("N_{Clusters}"); + fHistITSMultiplicity -> Sumw2(); + + fHistCentrality = new TH1D("fHistCentrality", Form("%s Centrality",fCentMethod.Data()), 300, -100, 200); + fHistCentrality -> SetXTitle("Centrality [%]"); + fHistCentrality -> Sumw2(); + + fOutputList -> Add(fHistV0Multiplicity); + fOutputList -> Add(fHistITSMultiplicity); + fOutputList -> Add(fHistCentrality); + + fHistEvStat = new TH1D("fHistEvStat","Event cuts statistics",20,-0.5,19.5); + fHistEvStat->SetXTitle("Cut index"); + fOutputList->Add(fHistEvStat); + + const Int_t kNZvtxBins = 10; + // bins for further buffers are shifted by 100 cm + Double_t vertexBins[kNZvtxBins+1] = { -10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10 }; + Int_t nZvtxBins = kNZvtxBins; + Double_t* zvtxbin = vertexBins; + + fPoolMgr = new AliEventPoolManager(1000, 20000, fNbinsCent, (Double_t*)fCentAxis->GetXbins()->GetArray(), nZvtxBins, zvtxbin); + + PostData(1, fOutputList); + +} + +//==================================================================================================================================================== + +void AliAnalysisTaskDiMuonCorrelations::UserExec(Option_t *) { + + fAOD = dynamic_cast(InputEvent()); + if (!fAOD) return; + + Int_t cutIndex = 0; + fHistEvStat->Fill(cutIndex++); + // Trigger selection + if (!(IsTriggerFired())) return; + fHistEvStat->Fill(cutIndex++); + + fHistV0Multiplicity -> Fill(GetV0Multiplicity()); + fHistITSMultiplicity -> Fill(GetITSMultiplicity()); + + Int_t centBin = GetCentBin(); + if (centBin<0) return; + fHistEvStat->Fill(cutIndex++); + Double_t percentile = fAOD->GetCentrality()->GetCentralityPercentile(fCentMethod.Data()); + fHistCentrality->Fill(percentile); + + // Vertex selection + const AliAODVertex* trkVtx = fAOD->GetPrimaryVertex(); + if (!trkVtx || trkVtx->GetNContributors()<=0) return; + TString vtxTtl = trkVtx->GetTitle(); + if (!vtxTtl.Contains("VertexerTracks")) return; + fHistEvStat->Fill(cutIndex++); + Double_t zvtx = trkVtx->GetZ(); + const AliAODVertex* spdVtx = fAOD->GetPrimaryVertexSPD(); + if (spdVtx->GetNContributors()<=0) return; + TString vtxTyp = spdVtx->GetTitle(); + Double_t cov[6]={0}; + spdVtx->GetCovarianceMatrix(cov); + Double_t zRes = TMath::Sqrt(cov[5]); + if (vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) return; + if (TMath::Abs(spdVtx->GetZ() - trkVtx->GetZ())>0.5) return; + fHistEvStat->Fill(cutIndex++); + + if (TMath::Abs(zvtx) > 10.) return; + fHistEvStat->Fill(cutIndex++); + + TObjArray *tracksMuonArm = GetAcceptedTracksMuonArm(fAOD); + if (tracksMuonArm->GetEntriesFast() < 2) { + delete tracksMuonArm; + return; + } + fHistEvStat->Fill(cutIndex++); + + fHistNMuons_vs_NMuons[centBin] -> Fill(tracksMuonArm->GetEntries(), tracksMuonArm->GetEntries()); + + AliDebug(1, Form("Single Event analysis : nTracksMuonArm = %4d", tracksMuonArm->GetEntries())); + + // Same event + for (Int_t iTrMuon1=0; iTrMuon1GetEntriesFast(); iTrMuon1++) { + fMuonTrack[0] = (AliAODTrack*) tracksMuonArm->At(iTrMuon1); + fHistSingleMuonsPt[centBin]->Fill(fMuonTrack[0]->Pt()); + fHistSingleMuonsEtaPt[centBin]->Fill(fMuonTrack[0]->Pt(),fMuonTrack[0]->Eta()); + for (Int_t iTrMuon2=iTrMuon1+1; iTrMuon2GetEntriesFast(); iTrMuon2++) { + fMuonTrack[1] = (AliAODTrack*) tracksMuonArm -> At(iTrMuon2); + FillHistograms(centBin, kSingleEvent); + } + } + + // Mixed event + AliEventPool* pool = fPoolMgr->GetEventPool(percentile, zvtx); + //pool->PrintInfo(); + if (pool->IsReady() || pool->NTracksInPool() > 2000 || pool->GetCurrentNEvents() >= 5) { + for (Int_t jMix=0; jMixGetCurrentNEvents(); jMix++) { + TObjArray *mixedTracks = pool->GetEvent(jMix); // Cvetan, here we would need to retrieve the MUON tracks of the event we mix... + fHistNMuons_vs_NMuons_Mixed[centBin]->Fill(mixedTracks->GetEntriesFast(), tracksMuonArm->GetEntriesFast()); + for (Int_t iTrMuon1=0; iTrMuon1GetEntriesFast(); iTrMuon1++) { + fMuonTrack[0] = (AliAODTrack*) tracksMuonArm->At(iTrMuon1); + fHistSingleMuonsPt_Mixed[centBin]->Fill(fMuonTrack[0]->Pt()); + fHistSingleMuonsEtaPt_Mixed[centBin]->Fill(fMuonTrack[0]->Pt(),fMuonTrack[0]->Eta()); + for (Int_t iTrMuon2=0; iTrMuon2GetEntriesFast(); iTrMuon2++) { + fMuonTrack[1] = (AliAODTrack*) mixedTracks -> At(iTrMuon2); + FillHistograms(centBin, kMixedEvent); + } + } + } + } + // pool->UpdatePool(tracksCentralBarrel); // Cvetan, I think I do not fully understand what this line does, and how it should be modified in a Di-Muon analysis + + delete tracksMuonArm; + + PostData(1, fOutputList); + +} + +//==================================================================================================================================================== + +void AliAnalysisTaskDiMuonCorrelations::FillHistograms(Int_t centrality, Int_t option) { + + Int_t ptBinTrackMuon1 = fPtAxis -> FindBin(fMuonTrack[0]->Pt()); + Int_t ptBinTrackMuon2 = fPtAxis -> FindBin(fMuonTrack[1]->Pt()); + + if (ptBinTrackMuon1<1 || ptBinTrackMuon1>fNbinsPt || ptBinTrackMuon2<1 || ptBinTrackMuon2>fNbinsPt) return; + + Double_t deltaPhi = fMuonTrack[0]->Phi() - fMuonTrack[1]->Phi(); + if (deltaPhi > 1.5*TMath::Pi()) deltaPhi -= TMath::TwoPi(); + if (deltaPhi < -0.5*TMath::Pi()) deltaPhi += TMath::TwoPi(); + + if (option==kSingleEvent) fHistDeltaPhi[centrality][ptBinTrackMuon1-1][ptBinTrackMuon2-1] -> Fill(TMath::RadToDeg()*deltaPhi); + else if (option==kMixedEvent) fHistDeltaPhiMix[centrality][ptBinTrackMuon1-1][ptBinTrackMuon2-1] -> Fill(TMath::RadToDeg()*deltaPhi); + + if (option==kSingleEvent) fHistEtaDeltaPhi[centrality][ptBinTrackMuon1-1][ptBinTrackMuon2-1] -> Fill(TMath::RadToDeg()*deltaPhi,fMuonTrack[0]->Eta()); + else if (option==kMixedEvent) fHistEtaDeltaPhiMix[centrality][ptBinTrackMuon1-1][ptBinTrackMuon2-1] -> Fill(TMath::RadToDeg()*deltaPhi,fMuonTrack[0]->Eta()); + + if (option==kSingleEvent) fHistTracksEtavsEta[centrality]->Fill(fMuonTrack[0]->Eta(),fMuonTrack[1]->Eta()); + else if (option==kMixedEvent) fHistTracksEtavsEta_Mixed[centrality]->Fill(fMuonTrack[0]->Eta(),fMuonTrack[1]->Eta()); + +} + +//==================================================================================================================================================== + +Bool_t AliAnalysisTaskDiMuonCorrelations::IsTriggerFired() { + + Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kINT7); + + return isSelected; +} + +//==================================================================================================================================================== + +Float_t AliAnalysisTaskDiMuonCorrelations::GetV0Multiplicity() { + + Float_t multiplicity=0; + for (Int_t iChannel=0; iChannel<64; iChannel++) multiplicity += fAOD->GetVZEROData()->GetMultiplicity(iChannel); + return multiplicity; + +} + +//==================================================================================================================================================== + +TObjArray* AliAnalysisTaskDiMuonCorrelations::GetAcceptedTracksMuonArm(AliAODEvent *aodEvent) { + + // fills the array of muon tracks that pass the cuts + + TObjArray *tracks = new TObjArray; + tracks->SetOwner(kFALSE); + + Int_t nTracks = aodEvent->GetNTracks(); + + AliAODTrack *track = 0; + + for (Int_t iTrack=0; iTrackGetTrack(iTrack); + if (track->IsMuonTrack() && track->GetMatchTrigger()>=fTriggerMatchLevelMuon) { + if (track->Chi2perNDF() < fMaxChi2Muon) { + if (track->Eta() > fMinEtaMuon && track->Eta() < fMaxEtaMuon) { + tracks->Add(new AliAODTrack(*track)); + } + } + } + } + + return tracks; + +} + +//==================================================================================================================================================== + +void AliAnalysisTaskDiMuonCorrelations::SetCentBinning(Int_t nBins, Double_t *limits) { + + if (nBins>fNMaxBinsCentrality) { + AliInfo(Form("WARNING : only %d centrality bins (out of the %d proposed) will be considered",fNMaxBinsCentrality,nBins)); + nBins = fNMaxBinsCentrality; + } + if (nBins<=0) { + AliInfo("WARNING : at least one centrality bin must be considered"); + nBins = 1; + } + + fNbinsCent = nBins; + fCentAxis = new TAxis(fNbinsCent, limits); + +} + +//==================================================================================================================================================== + +void AliAnalysisTaskDiMuonCorrelations::SetPtBinning(Int_t nBins, Double_t *limits) { + + if (nBins>fNMaxBinsPt) { + AliInfo(Form("WARNING : only %d pt bins (out of the %d proposed) will be considered",fNMaxBinsPt,nBins)); + nBins = fNMaxBinsPt; + } + if (nBins<=0) { + AliInfo("WARNING : at least one pt bin must be considered"); + nBins = 1; + } + + fNbinsPt = nBins; + fPtAxis = new TAxis(fNbinsPt, limits); + +} + +//==================================================================================================================================================== + +void AliAnalysisTaskDiMuonCorrelations::SetEtaBinning(Int_t nBins, Double_t *limits) { + + if (nBins>fNMaxBinsEta) { + AliInfo(Form("WARNING : only %d pt bins (out of the %d proposed) will be considered",fNMaxBinsEta,nBins)); + nBins = fNMaxBinsEta; + } + if (nBins<=0) { + AliInfo("WARNING : at least one pt bin must be considered"); + nBins = 1; + } + + fEtaAxis = new TAxis(nBins, limits); + +} + +//==================================================================================================================================================== + +Int_t AliAnalysisTaskDiMuonCorrelations::GetCentBin() { + + Double_t percentile = fAOD->GetCentrality()->GetCentralityPercentile(fCentMethod.Data()); + + Int_t bin = fCentAxis->FindBin(percentile) - 1; + if (bin >= fNbinsCent) bin = -1; + return bin; + +} + +//==================================================================================================================================================== + +Double_t AliAnalysisTaskDiMuonCorrelations::GetITSMultiplicity() { + + Double_t multiplicity = fAOD->GetHeader()->GetNumberOfITSClusters(1); + + return multiplicity; + +} + +//==================================================================================================================================================== + +void AliAnalysisTaskDiMuonCorrelations::Terminate(Option_t *) { + + fOutputList = dynamic_cast (GetOutputData(1)); + if (!fOutputList) { + printf("ERROR: Output list not available\n"); + return; + } + +} + +//==================================================================================================================================================== + diff --git a/PWGCF/Correlations/DPhi/MuonHadron/AliAnalysisTaskDiMuonCorrelations.h b/PWGCF/Correlations/DPhi/MuonHadron/AliAnalysisTaskDiMuonCorrelations.h new file mode 100644 index 00000000000..2b8bb3cfab2 --- /dev/null +++ b/PWGCF/Correlations/DPhi/MuonHadron/AliAnalysisTaskDiMuonCorrelations.h @@ -0,0 +1,100 @@ +#ifndef AliAnalysisTaskDiMuonCorrelations_H +#define AliAnalysisTaskDiMuonCorrelations_H + +#include "AliAnalysisTaskSE.h" +#include "TString.h" +#include "AliAODEvent.h" +#include "AliAODTrack.h" +#include "TAxis.h" +#include "TH1D.h" +#include "TClonesArray.h" +#include "TH2D.h" + +class AliEventPoolManager; + +//==================================================================================================================================================== + +class AliAnalysisTaskDiMuonCorrelations : public AliAnalysisTaskSE { + + public: + + enum {kSingleEvent, kMixedEvent}; + + AliAnalysisTaskDiMuonCorrelations(); + AliAnalysisTaskDiMuonCorrelations(const Char_t *name); + virtual ~AliAnalysisTaskDiMuonCorrelations(); + + virtual void UserCreateOutputObjects(); + virtual void UserExec(Option_t *option); + virtual void Terminate(Option_t *); + + // ------------- Cuts ----------------- + + void SetTriggerMatchLevelMuon(Short_t level) { fTriggerMatchLevelMuon = level; } + // void SetMaxChi2Muon(Double_t chi2Max) { fMaxChi2Muon = chi2Max; } + void SetEtaRangeMuon (Double_t etaMin, Double_t etaMax) { fMinEtaMuon = etaMin; fMaxEtaMuon = etaMax; } + + // ------------- Analysis ------------- + + Float_t GetV0Multiplicity(); + Double_t GetITSMultiplicity(); + Bool_t IsTriggerFired(); + TObjArray* GetAcceptedTracksMuonArm(AliAODEvent *aodEvent); + void SetPtBinning(Int_t nBins, Double_t *limits); + void SetCentBinning(Int_t nBins, Double_t *limits); + void SetEtaBinning(Int_t nBins, Double_t *limits); + void SetCentMethod(const Char_t *method) { fCentMethod = method; } + void FillHistograms(Int_t centrality, Int_t option); + Int_t GetCentBin(); + + private: + + static const Int_t fNMaxBinsCentrality = 20; + static const Int_t fNMaxBinsPt = 10; + static const Int_t fNMaxBinsEta = 10; + + AliAODEvent *fAOD; //! + AliEventPoolManager *fPoolMgr; //! event pool manager + AliAODTrack *fMuonTrack[2]; //! + + Double_t fMaxChi2Muon, fMinEtaMuon, fMaxEtaMuon; + Short_t fTriggerMatchLevelMuon; + + Int_t fNbinsCent, fNbinsPt; + + TAxis *fCentAxis; + TAxis *fPtAxis; + TAxis *fEtaAxis; + + TH1D *fHistDeltaPhi[fNMaxBinsCentrality][fNMaxBinsPt][fNMaxBinsPt]; //! + TH1D *fHistDeltaPhiMix[fNMaxBinsCentrality][fNMaxBinsPt][fNMaxBinsPt]; //! + TH2D *fHistEtaDeltaPhi[fNMaxBinsCentrality][fNMaxBinsPt][fNMaxBinsPt]; //! + TH2D *fHistEtaDeltaPhiMix[fNMaxBinsCentrality][fNMaxBinsPt][fNMaxBinsPt]; //! + TH2D *fHistNMuons_vs_NMuons[fNMaxBinsCentrality]; //! + TH2D *fHistNMuons_vs_NMuons_Mixed[fNMaxBinsCentrality]; //! + TH2D *fHistTracksEtavsEta[fNMaxBinsCentrality]; //! + TH2D *fHistTracksEtavsEta_Mixed[fNMaxBinsCentrality]; //! + TH1D *fHistSingleMuonsPt[fNMaxBinsCentrality]; //! + TH1D *fHistSingleMuonsPt_Mixed[fNMaxBinsCentrality]; //! + TH2D *fHistSingleMuonsEtaPt[fNMaxBinsCentrality]; //! + TH2D *fHistSingleMuonsEtaPt_Mixed[fNMaxBinsCentrality]; //! + + TH1D *fHistV0Multiplicity; //! + TH1D *fHistITSMultiplicity; //! + TH1D *fHistCentrality; //! + TH1D *fHistEvStat; //! + + TString fCentMethod; + + TList *fOutputList; //! + + AliAnalysisTaskDiMuonCorrelations(const AliAnalysisTaskDiMuonCorrelations&);//not implimented + AliAnalysisTaskDiMuonCorrelations& operator=(const AliAnalysisTaskDiMuonCorrelations&);//not implimnted + + ClassDef(AliAnalysisTaskDiMuonCorrelations, 2) // example of analysis + +}; + +//==================================================================================================================================================== + +#endif diff --git a/PWGCF/PWGCFCorrelationsDPhiLinkDef.h b/PWGCF/PWGCFCorrelationsDPhiLinkDef.h index 71d0e037bcb..446a8274242 100644 --- a/PWGCF/PWGCFCorrelationsDPhiLinkDef.h +++ b/PWGCF/PWGCFCorrelationsDPhiLinkDef.h @@ -23,4 +23,5 @@ #pragma link C++ class AliAnalysisTaskLongRangeCorrelations+; #pragma link C++ class LRCParticle; #pragma link C++ class AliAnalysisTaskMuonHadronCorrelations+; +#pragma link C++ class AliAnalysisTaskDiMuonCorrelations+; #endif -- 2.43.0