1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
17 Simple use case for mixed event analysis
19 Delta_phi correlation analysis is performed on charged tracks
20 for same and mixed events
21 Author: andreas.morsch@cern.ch
28 #include "AliAODEvent.h"
29 #include "AliAODTrack.h"
30 #include "AliAnalysisTaskPhiCorr.h"
31 #include "AliMultiEventInputHandler.h"
33 ClassImp(AliAnalysisTaskPhiCorr)
35 //________________________________________________________________________
36 AliAnalysisTaskPhiCorr::AliAnalysisTaskPhiCorr(const char *name)
37 : AliAnalysisTaskME(name), fHists(0), fHistDphiCO(0), fHistDphiUC(0)
40 // Define input and output slots here
41 DefineOutput(1, TList::Class());
45 //________________________________________________________________________
46 void AliAnalysisTaskPhiCorr::UserCreateOutputObjects()
51 fHistDphiCO = new TH1F("fHistDPhiCO", "#Delta #Phi distribution", 64, 0., 3.2);
52 fHistDphiCO->GetXaxis()->SetTitle("#Delta#Phi [rad]");
53 fHistDphiCO->GetYaxis()->SetTitle("dN/d#Phi");
54 fHistDphiCO->SetMarkerStyle(kFullCircle);
56 fHistDphiUC = new TH1F("fHistDPhiUC", "#Delta #Phi distribution", 64, 0., 3.2);
57 fHistDphiUC->GetXaxis()->SetTitle("#Delta#Phi [rad]");
58 fHistDphiUC->GetYaxis()->SetTitle("dN/d#Phi");
59 fHistDphiUC->SetMarkerStyle(kOpenCircle);
62 fHists->Add(fHistDphiCO);
63 fHists->Add(fHistDphiUC);
66 //________________________________________________________________________
67 void AliAnalysisTaskPhiCorr::UserExec(Option_t *)
70 // Called for each event
71 // Uncorrelated tracks
72 Int_t nev = fInputHandler->GetBufferSize();
73 Float_t wgt = 1./(nev*(nev-1));
77 for (Int_t iev = 0; iev < nev; iev++) {
78 AliAODEvent* aod = (AliAODEvent*) GetEvent(iev);
79 fMixedEvent.AddEvent(aod);
82 Int_t ntrack = fMixedEvent.GetNumberOfTracks();
84 for (Int_t itr = 0; itr < ntrack -1; itr++) {
85 for (Int_t jtr = itr+1; jtr < ntrack; jtr++) {
86 AliVParticle* track1 = fMixedEvent.GetTrack(itr);
87 AliVParticle* track2 = fMixedEvent.GetTrack(jtr);
88 Int_t iev1 = fMixedEvent.EventIndex(itr);
89 Int_t iev2 = fMixedEvent.EventIndex(jtr);
91 Float_t phi1 = track1->Phi();
92 Float_t phi2 = track2->Phi();
93 Float_t dphi = TMath::Abs(phi1 - phi2);
94 if (dphi > TMath::Pi()) dphi = 2. * TMath::Pi() - dphi;
96 fHistDphiUC->Fill(dphi, wgt);
98 fHistDphiCO->Fill(dphi, 0.5);
107 //________________________________________________________________________
108 void AliAnalysisTaskPhiCorr::Terminate(Option_t *)
110 // Draw result to the screen
111 // Called once at the end of the query
113 TCanvas *c1 = new TCanvas("AliAnalysisTaskPt","Pt",10,10,510,510);
114 c1->cd(1)->SetLogy();
115 fHistDphiCO->DrawCopy("E");
116 fHistDphiUC->DrawCopy("Esame");