6 #include "AliAODEvent.h"
7 #include "AliAODTrack.h"
11 #include "AliAnalysisTaskMuonHadronCorrelations.h"
13 ClassImp(AliAnalysisTaskMuonHadronCorrelations)
15 //====================================================================================================================================================
17 AliAnalysisTaskMuonHadronCorrelations::AliAnalysisTaskMuonHadronCorrelations() :
20 fTracksCentralBarrel(0x0),
24 fFilterBitCentralBarrel(0),
25 fMaxEtaCentralBarrel(1.0),
26 fMaxChi2Muon(9999999999.),
28 fMaxRAbsMuon(9999999999.),
29 fTriggerMatchLevelMuon(0),
31 fIsTriggerSet(kFALSE),
37 fHistV0Multiplicity(0x0),
38 fHistITSMultiplicity(0x0),
45 // Default constructor
47 for (Int_t iCent=0; iCent<fNMaxBinsCentrality; iCent++) {
48 for (Int_t iPtBinCB=0; iPtBinCB<fNMaxBinsPt; iPtBinCB++) {
49 for (Int_t iPtBinMA=0; iPtBinMA<fNMaxBinsPt; iPtBinMA++) {
50 fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA] = NULL;
51 fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA] = NULL;
54 fHistNTracksCB_vs_NTracksMA[iCent] = NULL;
60 //====================================================================================================================================================
62 AliAnalysisTaskMuonHadronCorrelations::AliAnalysisTaskMuonHadronCorrelations(const char *name) :
63 AliAnalysisTaskSE(name),
65 fTracksCentralBarrel(0x0),
69 fFilterBitCentralBarrel(0),
70 fMaxEtaCentralBarrel(1.0),
71 fMaxChi2Muon(9999999999.),
73 fMaxRAbsMuon(9999999999.),
74 fTriggerMatchLevelMuon(0),
76 fIsTriggerSet(kFALSE),
82 fHistV0Multiplicity(0x0),
83 fHistITSMultiplicity(0x0),
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;
99 fHistNTracksCB_vs_NTracksMA[iCent] = NULL;
102 // Define input and output slots here
103 DefineOutput(1, TList::Class());
107 //====================================================================================================================================================
109 AliAnalysisTaskMuonHadronCorrelations::~AliAnalysisTaskMuonHadronCorrelations() {
118 for (Int_t iCent=0; iCent<fNMaxBinsCentrality; iCent++) {
119 for (Int_t iPtBinCB=0; iPtBinCB<fNMaxBinsPt; iPtBinCB++) {
120 for (Int_t iPtBinMA=0; iPtBinMA<fNMaxBinsPt; iPtBinMA++) {
121 delete fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA];
122 delete fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA];
125 delete fHistNTracksCB_vs_NTracksMA[iCent];
130 //====================================================================================================================================================
132 void AliAnalysisTaskMuonHadronCorrelations::UserCreateOutputObjects() {
134 fOutputList = new TList();
135 fOutputList->SetOwner(kTRUE);
137 for (Int_t iCent=0; iCent<fNbinsCent; iCent++) {
138 for (Int_t iPtBinCB=0; iPtBinCB<fNbinsPt; iPtBinCB++) {
139 for (Int_t iPtBinMA=0; iPtBinMA<fNbinsPt; iPtBinMA++) {
141 fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA] = new TH1D(Form("fHistDeltaPhi_Cent%02d_PtBin%02d_%02d",iCent,iPtBinCB,iPtBinMA),
142 Form("%d-%d %%, %3.1f<p_{T}^{TPC}<%3.1f, %3.1f<p_{T}^{Muon}<%3.1f",
143 Int_t(fCentAxis->GetBinLowEdge(iCent+1)),
144 Int_t(fCentAxis->GetBinUpEdge(iCent+1)),
145 fPtAxis->GetBinLowEdge(iPtBinCB+1),
146 fPtAxis->GetBinUpEdge(iPtBinCB+1),
147 fPtAxis->GetBinLowEdge(iPtBinMA+1),
148 fPtAxis->GetBinUpEdge(iPtBinMA+1)),
149 100, -0.5*TMath::RadToDeg()*TMath::Pi(), 1.5*TMath::RadToDeg()*TMath::Pi());
151 fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA] = new TH1D(Form("fHistDeltaPhiMix_Cent%02d_PtBin%02d_%02d",iCent,iPtBinCB,iPtBinMA),
152 Form("%d-%d %%, %3.1f<p_{T}^{TPC}<%3.1f, %3.1f<p_{T}^{Muon}<%3.1f MIXED",
153 Int_t(fCentAxis->GetBinLowEdge(iCent+1)),
154 Int_t(fCentAxis->GetBinUpEdge(iCent+1)),
155 fPtAxis->GetBinLowEdge(iPtBinCB+1),
156 fPtAxis->GetBinUpEdge(iPtBinCB+1),
157 fPtAxis->GetBinLowEdge(iPtBinMA+1),
158 fPtAxis->GetBinUpEdge(iPtBinMA+1)),
159 100, -0.5*TMath::RadToDeg()*TMath::Pi(), 1.5*TMath::RadToDeg()*TMath::Pi());
161 fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA] -> SetXTitle("#Delta#varphi [degrees]");
162 fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA] -> SetXTitle("#Delta#varphi [degrees]");
164 fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA] -> Sumw2();
165 fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA] -> Sumw2();
167 fOutputList -> Add(fHistDeltaPhi[iCent][iPtBinCB][iPtBinMA]);
168 fOutputList -> Add(fHistDeltaPhiMix[iCent][iPtBinCB][iPtBinMA]);
173 fHistNTracksCB_vs_NTracksMA[iCent] = new TH2D(Form("fHistNTracksCB_vs_NTracksMA_Cent%02d",iCent),
174 Form("%d-%d %%",Int_t(fCentAxis->GetBinLowEdge(iCent+1)),Int_t(fCentAxis->GetBinUpEdge(iCent+1))),
175 100, 0, 1000, 100, 0, 100);
176 fHistNTracksCB_vs_NTracksMA[iCent] -> SetXTitle("N_{tracks} Central Barrel");
177 fHistNTracksCB_vs_NTracksMA[iCent] -> SetYTitle("N_{tracks} Muon Arm");
179 fHistNTracksCB_vs_NTracksMA[iCent] -> Sumw2();
181 fOutputList -> Add(fHistNTracksCB_vs_NTracksMA[iCent]);
185 fHistV0Multiplicity = new TH1D("fHistV0Multiplicity", "V0 Multiplicity", 500, 0, 1000);
186 fHistV0Multiplicity -> SetXTitle("Multiplicity");
187 fHistV0Multiplicity -> Sumw2();
189 fHistITSMultiplicity = new TH1D("fHistITSMultiplicity", "ITS Multiplicity", 500, 0, 500);
190 fHistITSMultiplicity -> SetXTitle("N_{Clusters}");
191 fHistITSMultiplicity -> Sumw2();
193 fHistCentrality = new TH1D("fHistCentrality", Form("%s Centrality",fCentMethod.Data()), 300, -100, 200);
194 fHistCentrality -> SetXTitle("Centrality [%]");
195 fHistCentrality -> Sumw2();
197 fOutputList -> Add(fHistV0Multiplicity);
198 fOutputList -> Add(fHistITSMultiplicity);
199 fOutputList -> Add(fHistCentrality);
201 PostData(1, fOutputList);
205 //====================================================================================================================================================
207 void AliAnalysisTaskMuonHadronCorrelations::UserExec(Option_t *) {
209 AliDebug(2, Form("Single Event analysis : event %05d",fEvt));
211 fAOD = dynamic_cast<AliAODEvent *>(InputEvent());
215 if (!(IsTriggerFired())) return;
217 fHistV0Multiplicity -> Fill(GetV0Multiplicity());
218 fHistITSMultiplicity -> Fill(GetITSMultiplicity());
220 Int_t centBin = GetCentBin();
221 if (centBin<0) return;
223 fTracksCentralBarrel = GetAcceptedTracksCentralBarrel(fAOD);
224 fTracksMuonArm = GetAcceptedTracksMuonArm(fAOD);
226 fHistNTracksCB_vs_NTracksMA[centBin] -> Fill(fTracksCentralBarrel->GetEntries(), fTracksMuonArm->GetEntries());
228 AliDebug(1, Form("Single Event analysis : event %05d, nTracksCB = %4d, nTracksMA = %4d",fEvt, fTracksCentralBarrel->GetEntries(), fTracksMuonArm->GetEntries()));
230 for (Int_t iTrCB=0; iTrCB<fTracksCentralBarrel->GetEntriesFast(); iTrCB++) {
231 for (Int_t iTrMA=0; iTrMA<fTracksMuonArm->GetEntriesFast(); iTrMA++) {
232 fTrackCB = (AliAODTrack*) fTracksCentralBarrel -> At(iTrCB);
233 fTrackMA = (AliAODTrack*) fTracksMuonArm -> At(iTrMA);
234 FillHistograms(centBin, kSingleEvent);
238 delete fTracksCentralBarrel;
239 delete fTracksMuonArm;
241 PostData(1, fOutputList);
245 //====================================================================================================================================================
247 void AliAnalysisTaskMuonHadronCorrelations::FillHistograms(Int_t centrality, Int_t option) {
249 Int_t ptBinTrackCB = fPtAxis -> FindBin(fTrackCB->Pt());
250 Int_t ptBinTrackMA = fPtAxis -> FindBin(fTrackMA->Pt());
252 if (ptBinTrackCB<1 || ptBinTrackCB>fNbinsPt || ptBinTrackMA<1 || ptBinTrackMA>fNbinsPt) return;
254 Double_t deltaPhi = fTrackCB->Phi() - fTrackMA->Phi();
255 if (deltaPhi > 1.5*TMath::Pi()) deltaPhi -= TMath::TwoPi();
256 if (deltaPhi < -0.5*TMath::Pi()) deltaPhi += TMath::TwoPi();
258 if (option==kSingleEvent) fHistDeltaPhi[centrality][ptBinTrackCB-1][ptBinTrackMA-1] -> Fill(TMath::RadToDeg()*deltaPhi);
259 else if (option==kMixedEvent) fHistDeltaPhiMix[centrality][ptBinTrackCB-1][ptBinTrackMA-1] -> Fill(TMath::RadToDeg()*deltaPhi);
263 //====================================================================================================================================================
265 Bool_t AliAnalysisTaskMuonHadronCorrelations::IsTriggerFired() {
267 Bool_t v0and = ((fAOD->GetVZEROData()->GetV0ADecision()==1) && (fAOD->GetVZEROData()->GetV0CDecision()==1));
268 TString trigStr(fAOD->GetHeader()->GetFiredTriggerClasses());
269 return (trigStr.Contains(fTriggerWord) && v0and);
273 //====================================================================================================================================================
275 Float_t AliAnalysisTaskMuonHadronCorrelations::GetV0Multiplicity() {
277 Float_t multiplicity=0;
278 for (Int_t iChannel=0; iChannel<64; iChannel++) multiplicity += fAOD->GetVZEROData()->GetMultiplicity(iChannel);
283 //====================================================================================================================================================
285 TClonesArray* AliAnalysisTaskMuonHadronCorrelations::GetAcceptedTracksCentralBarrel(AliAODEvent *aodEvent) {
287 // fills the array of central barrel tracks that pass the cuts
289 TClonesArray *tracks = new TClonesArray("AliAODTrack");
291 Int_t nTracks = aodEvent->GetNTracks();
293 AliAODTrack *track = 0;
295 for (Int_t iTrack=0; iTrack<nTracks; iTrack++) {
296 track = aodEvent->GetTrack(iTrack);
297 if (track->TestFilterBit(fFilterBitCentralBarrel) && TMath::Abs(track->Eta())<fMaxEtaCentralBarrel) {
298 new ((*tracks)[tracks->GetEntries()]) AliAODTrack(*track);
306 //====================================================================================================================================================
308 TClonesArray* AliAnalysisTaskMuonHadronCorrelations::GetAcceptedTracksMuonArm(AliAODEvent *aodEvent) {
310 // fills the array of muon tracks that pass the cuts
312 TClonesArray *tracks = new TClonesArray("AliAODTrack");
314 Int_t nTracks = aodEvent->GetNTracks();
316 AliAODTrack *track = 0;
318 for (Int_t iTrack=0; iTrack<nTracks; iTrack++) {
319 track = aodEvent->GetTrack(iTrack);
320 if (track->IsMuonTrack() && track->GetMatchTrigger()>=fTriggerMatchLevelMuon) {
321 new ((*tracks)[tracks->GetEntries()]) AliAODTrack(*track);
329 //====================================================================================================================================================
331 void AliAnalysisTaskMuonHadronCorrelations::SetCentBinning(Int_t nBins, Double_t *limits) {
333 if (nBins>fNMaxBinsCentrality) {
334 AliInfo(Form("WARNING : only %d centrality bins (out of the %d proposed) will be considered",fNMaxBinsCentrality,nBins));
335 nBins = fNMaxBinsCentrality;
338 AliInfo("WARNING : at least one centrality bin must be considered");
343 fCentAxis = new TAxis(fNbinsCent, limits);
347 //====================================================================================================================================================
349 void AliAnalysisTaskMuonHadronCorrelations::SetMultBinning(Int_t nBins, Double_t *limits) {
351 if (nBins>fNMaxBinsCentrality) {
352 AliInfo(Form("WARNING : only %d centrality bins (out of the %d proposed) will be considered",fNMaxBinsCentrality,nBins));
353 nBins = fNMaxBinsCentrality;
356 AliInfo("WARNING : at least one centrality bin must be considered");
361 fMultAxis = new TAxis(fNbinsCent, limits);
365 //====================================================================================================================================================
367 void AliAnalysisTaskMuonHadronCorrelations::SetPtBinning(Int_t nBins, Double_t *limits) {
369 if (nBins>fNMaxBinsPt) {
370 AliInfo(Form("WARNING : only %d pt bins (out of the %d proposed) will be considered",fNMaxBinsPt,nBins));
374 AliInfo("WARNING : at least one pt bin must be considered");
379 fPtAxis = new TAxis(fNbinsPt, limits);
383 //====================================================================================================================================================
385 Int_t AliAnalysisTaskMuonHadronCorrelations::GetCentBin() {
387 if (fCentMethod.CompareTo("V0M")==0) return GetMultBin();
388 if (fCentMethod.CompareTo("CL1")==0) return GetMultBin();
394 //====================================================================================================================================================
396 Int_t AliAnalysisTaskMuonHadronCorrelations::GetMultBin() {
398 Double_t multiplicity = 0;
400 if (fCentMethod.CompareTo("V0M")==0) multiplicity = GetV0Multiplicity();
401 if (fCentMethod.CompareTo("CL1")==0) multiplicity = GetITSMultiplicity();
403 Int_t multBin = fMultAxis->FindBin(-1.*multiplicity);
404 if (0<multBin && multBin<=fNbinsCent) return multBin-1;
410 //====================================================================================================================================================
412 Double_t AliAnalysisTaskMuonHadronCorrelations::GetITSMultiplicity() {
414 Double_t multiplicity = 0.5 * (fAOD->GetHeader()->GetNumberOfITSClusters(0) + fAOD->GetHeader()->GetNumberOfITSClusters(1));
420 //====================================================================================================================================================
422 void AliAnalysisTaskMuonHadronCorrelations::Terminate(Option_t *) {
424 fOutputList = dynamic_cast<TList*> (GetOutputData(1));
426 printf("ERROR: Output list not available\n");
432 //====================================================================================================================================================