Adding analysis for Flat ESD case
[u/mrichter/AliRoot.git] / TPC / Calib / AliAnalysisTaskPt.cxx
CommitLineData
4d3da966 1#include "TChain.h"
2#include "TTree.h"
3#include "TH1F.h"
4#include "TCanvas.h"
5#include "TObjArray.h"
6
7#include "AliAnalysisTask.h"
8#include "AliAnalysisManager.h"
9
10#include "AliESDEvent.h"
11#include "AliVVevent.h"
12#include "AliVVtrack.h"
13#include "AliESDtrackCuts.h"
14#include "AliVEventHandler.h"
15#include "AliTPCseed.h"
16#include "AliTPCclusterMI.h"
17#include "AliVVfriendEvent.h"
18#include "AliVVfriendTrack.h"
19
20#include "AliAnalysisTaskPt.h"
21
22// example of an analysis task creating a p_t spectrum
23// Authors: Panos Cristakoglou, Jan Fiete Grosse-Oetringhaus, Christian Klein-Boesing
24
25ClassImp(AliAnalysisTaskPt)
26
27//________________________________________________________________________
28AliAnalysisTaskPt::AliAnalysisTaskPt(const char *name)
828c4c08 29: AliAnalysisTask(name, ""), fESD(0), fESDfriend(0), fHistPt(0), fCuts(0), fEv(0), fHistQ(0), fListOut(0), fUseFriends(kFALSE), fHistNTPCCl(0), fHistNESDtracks(0), fHistNESDfriendtracks(0)
30
4d3da966 31{
32 // Constructor
33
34 // Define input and output slots here
35 // Input slot #0 works with a TChain
36 DefineInput(0, TChain::Class());
37 // Output slot #0 writes into a TH1 container
38 DefineOutput(0, TList::Class());
39}
40
41//________________________________________________________________________
42void AliAnalysisTaskPt::ConnectInputData(Option_t *)
43{
44 // Connect ESD or AOD here
45 // Called once
46
47 printf("----> AliAnalysisTaskPt::ConnectInputData\n");
48 TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
49 if (!tree) {
50 Printf("ERROR: Could not read chain from input slot 0");
51 } else {
52 // Disable all branches and enable only the needed ones
53 // The next two lines are different when data produced as AliESDEvent is read
54 /*
55 tree->SetBranchStatus("*", kFALSE);
56 tree->SetBranchStatus("fTracks.*", kTRUE);
57 */
58
59 AliVEventHandler *esdH = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
60
61 if (!esdH) {
62 Printf("ERROR: Could not get ESDInputHandler");
63 } else {
64 Printf("----> AliAnalysisTaskPt::ConnectInputData Getting the event from handler %p", esdH);
65 //fESD = dynamic_cast<AliESDEvent*>(esdH->GetEvent());
a453390d 66 fESD = esdH->GetVVEvent();
67 if (fUseFriends){
68 fESDfriend = esdH->GetVVFriendEvent();
69 }
4d3da966 70 }
71 if (!fESD) {
72 Printf("ERROR, no ESD event");
73 }
a453390d 74 if (fUseFriends && !fESDfriend) {
4d3da966 75 Printf("ERROR, no ESD friend");
76 }
77 }
78
79 Printf("fESD = %p, fESDfriend = %p", fESD, fESDfriend);
80 printf("<---- AliAnalysisTaskPt::ConnectInputData\n");
81}
82
83//________________________________________________________________________
84void AliAnalysisTaskPt::CreateOutputObjects()
85{
86 // Create histograms
87 // Called once
88
89 fListOut = new TList();
90 fListOut->SetOwner();
91 fListOut->SetName("listHistos");
92
93 fHistPt = new TH1F("fHistPt", "P_{T} distribution", 15, 0.1, 3.1);
94 fHistPt->GetXaxis()->SetTitle("P_{T} (GeV/c)");
95 fHistPt->GetYaxis()->SetTitle("dN/dP_{T} (c/GeV)");
96 fHistPt->SetMarkerStyle(kFullCircle);
97
98 fHistQ = new TH1F("fHistQ", "TPC clusters: Q distribution", 1000, 0, 10000);
99 fHistQ->GetXaxis()->SetTitle("Q");
100 fHistQ->GetYaxis()->SetTitle("dN/dQ");
101 fHistQ->SetMarkerStyle(kFullCircle);
102
828c4c08 103 fHistNTPCCl = new TH1F("fHistNTPCCl", "Number of TPC clusters", 160, -0.5, 159.5);
104 fHistNTPCCl->GetXaxis()->SetTitle("n. TPC Cl.");
105 fHistNTPCCl->GetYaxis()->SetTitle("dN/d(n. TPC Cl)");
106 fHistNTPCCl->SetMarkerStyle(kFullCircle);
107
108 fHistNESDtracks = new TH1F("fHistNESDtracks", "Number of ESD friend tracks", 1000, -0.5, 999.5);
109 fHistNESDtracks->GetXaxis()->SetTitle("n. ESD friend tracks");
110 fHistNESDtracks->GetYaxis()->SetTitle("dN/d(n. ESD friend tracks)");
111 fHistNESDtracks->SetMarkerStyle(kFullCircle);
112
113 fHistNESDfriendtracks = new TH1F("fHistNESDfriendtracks", "Number of ESD tracks", 1000, -0.5, 999.5);
114 fHistNESDfriendtracks->GetXaxis()->SetTitle("n. ESD tracks");
115 fHistNESDfriendtracks->GetYaxis()->SetTitle("dN/d(n. ESD tracks)");
116 fHistNESDfriendtracks->SetMarkerStyle(kFullCircle);
117
4d3da966 118 fListOut->Add(fHistPt);
119 fListOut->Add(fHistQ);
828c4c08 120 fListOut->Add(fHistNTPCCl);
121 fListOut->Add(fHistNESDtracks);
122 fListOut->Add(fHistNESDfriendtracks);
4d3da966 123
124 PostData(0, fListOut);
125
126 fCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(1);
127}
128
129//________________________________________________________________________
130void AliAnalysisTaskPt::Exec(Option_t *)
131{
132 // Main loop
133 // Called for each event
134
135 if (!fESD) {
136 Printf("ERROR: fESD not available");
137 return;
138 }
139 if (!fESDfriend) {
140 Printf("ERROR: fESDfriend not available");
828c4c08 141 if (fUseFriends){
142 return;
143 }
4d3da966 144 }
145
828c4c08 146 Int_t nESDtracks = fESD->GetNumberOfTracks();
147 Int_t nESDfriendtracks = 0;
148 if (fUseFriends) nESDfriendtracks = fESDfriend->GetNumberOfTracks();
149 Printf("There are %d tracks in this event", nESDtracks);
150 Printf("... and there are %d friends in this event", nESDfriendtracks);
151
152 fHistNESDtracks->Fill(nESDtracks);
153 fHistNESDfriendtracks->Fill(nESDfriendtracks);
4d3da966 154
155 // Track loop to fill a pT spectrum
828c4c08 156 for (Int_t iTracks = 0; iTracks < nESDtracks; iTracks++) {
157 Printf("Checking track %d: Note that with Flat, the GetTrack is not yet implemented!!!", iTracks);
6a33e0e9 158 AliVVtrack* track = fESD->GetVVTrack(iTracks);
4d3da966 159 if (!track) {
160 Printf("ERROR: Could not receive track %d", iTracks);
161 continue;
162 }
6a33e0e9 163 Printf("track %d has pt = %f", iTracks, track->GetPt());
164 fHistPt->Fill(track->GetPt());
828c4c08 165 fHistNTPCCl->Fill(track->GetTPCNcls());
4d3da966 166 } //track loop
167
168
a453390d 169 if (fUseFriends){
170 // Friend Track loop
828c4c08 171 for (Int_t iFriend = 0; iFriend < nESDfriendtracks; iFriend++) {
a453390d 172 //Printf("Getting friend %d", iFriend);
173 const AliVVfriendTrack* friendTrack = fESDfriend->GetTrack(iFriend);
174 if (!friendTrack) {
175 Printf("ERROR: Could not receive track %d", iFriend);
176 continue;
6a33e0e9 177 }
178
179 AliTPCseed seed;
180 Int_t err = friendTrack->GetTPCseed( seed );
181 if( err==0 ){
182 Printf("Found TPC seed" );
183 for (Int_t irow = 0; irow < 160; irow++){
184 AliTPCclusterMI* cluMI = seed.GetClusterPointer(irow);
185 if (cluMI){
186 Printf("Found cluster at row %d", irow);
187 Float_t q = cluMI->GetQ();
188 Printf("Q = %f", q);
189 fHistQ->Fill(q);
a453390d 190 }
6a33e0e9 191 else {
192 Printf("Row %d does not contain clusters", irow);
193 }
a453390d 194 }
6a33e0e9 195 }
196 else {
197 //Printf("Schade... seed is %p", seed);
4d3da966 198 }
a453390d 199 }
4d3da966 200 }
201
202 // Post output data.
203 PostData(0, fListOut);
204 fEv++;
205}
206
207//________________________________________________________________________
208void AliAnalysisTaskPt::Terminate(Option_t *)
209{
210 // Draw result to the screen
211 // Called once at the end of the query
212
213 Printf("Terminate called: fESD = %p", fESD);
214
215 fListOut = dynamic_cast<TList*> (GetOutputData(0));
216
217 if (fListOut) {
218 fHistPt = dynamic_cast<TH1F*>(fListOut->FindObject("fHistPt"));
219 if (!fHistPt) {
220 Printf("ERROR: fHistPt not available");
221 return;
222 }
223
224 TCanvas *c1 = new TCanvas("AliAnalysisTaskPt","Pt",10,10,510,510);
225 c1->cd(1)->SetLogy();
226 fHistPt->DrawCopy("E");
227 }
228 else {
229 Printf("In Terminate: no TList found");
230 }
231
232}