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