uuhhhhaaa - what did I change? Added the run number and fixed coding conventions ;-)
[u/mrichter/AliRoot.git] / PWG0 / TPC / AliROCESDAnalysisSelector.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
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  **************************************************************************/
15
16 /* $Id$ */
17
18 // The ESD is available as member fESD
19 //
20 // The Process function is nearly empty. Implement your analysis there and look at the other listed below functions you
21 // might need.
22 //
23 // The following methods can be overrriden. Please do not forgot to call the base class function.
24 //
25 //    Begin():        called everytime a loop on the tree starts,
26 //                    a convenient place to create your histograms.
27 //    SlaveBegin():   called after Begin(), when on PROOF called only on the
28 //                    slave servers.
29 //    Init():         called for each new tree. Enable/Disable branches here.
30 //    Process():      called for each event, in this function you decide what
31 //                    to read and fill your histograms.
32 //    SlaveTerminate: called at the end of the loop on the tree, when on PROOF
33 //                    called only on the slave servers.
34 //    Terminate():    called at the end of the loop on the tree,
35 //                    a convenient place to draw/fit your histograms.
36 //
37 //  Author: Jan.Fiete.Grosse-Oetringhaus@cern.ch
38
39 #include "AliROCESDAnalysisSelector.h"
40
41 #include <AliLog.h>
42 #include <AliESD.h>
43 #include <AliESDfriend.h>
44 #include <../TPC/AliTPCclusterMI.h>
45 #include <../TPC/AliTPCseed.h>
46
47 #include <TFile.h>
48 #include <TTree.h>
49 #include <TCanvas.h>
50
51 #include "TPC/AliTPCClusterHistograms.h"
52
53 ClassImp(AliROCESDAnalysisSelector)
54
55 AliROCESDAnalysisSelector::AliROCESDAnalysisSelector() :
56   AliSelector(),
57   fESDfriend(0)
58 {
59   //
60   // Constructor. Initialization of pointers
61   //
62   
63   for (Int_t i=0; i<kTPCSectors; i++)
64     fClusterHistograms[i] = 0;
65 }
66
67 AliROCESDAnalysisSelector::~AliROCESDAnalysisSelector()
68 {
69   //
70   // Destructor
71   //
72 }
73
74 void AliROCESDAnalysisSelector::SlaveBegin(TTree* tree)
75 {
76   //
77   
78   AliSelector::SlaveBegin(tree);
79
80
81 void AliROCESDAnalysisSelector::Init(TTree *tree)
82 {
83   // The Init() function is called when the selector needs to initialize
84   // a new tree or chain. Typically here the branch addresses of the tree
85   // will be set. It is normaly not necessary to make changes to the
86   // generated code, but the routine can be extended by the user if needed.
87   // Init() will be called many times when running with PROOF.
88
89   AliSelector::Init(tree);
90
91   // Set branch address
92   if (tree) {
93     tree->SetBranchStatus("*", 0);
94     tree->SetBranchStatus("fTracks.*", 1);
95     tree->SetBranchStatus("fTimeStamp", 1);
96     tree->SetBranchAddress("ESDfriend", &fESDfriend);
97   }
98   if (fESDfriend != 0)
99     AliDebug(AliLog::kInfo, "INFO: Found ESDfriend branch in chain.");
100 }
101
102 Bool_t AliROCESDAnalysisSelector::Process(Long64_t entry)
103 {
104   //
105   // Implement your analysis here. Do not forget to call the parent class Process by
106   // if (AliSelector::Process(entry) == kFALSE)
107   //   return kFALSE;
108   //
109
110   //  AliDebug(AliLog::kInfo, Form("Processing event %d \n", entry));
111
112
113   if (AliSelector::Process(entry) == kFALSE)
114     return kFALSE;
115
116   // Check prerequisites
117   if (!fESD)
118   {
119     AliDebug(AliLog::kError, "ESD branch not available");
120     return kFALSE;
121   }
122
123   // Check prerequisites
124   if (!fESDfriend)
125   {
126     AliDebug(AliLog::kError, "ESDfriend branch not available");
127     return kFALSE;
128   }
129
130   //  printf(Form(" event number: %d  ... time stamp: %d \n", fESD->GetEventNumber(), fESD->GetTimeStamp()));
131
132   fESD->SetESDfriend(fESDfriend);
133
134   Int_t nTracks = fESD->GetNumberOfTracks();
135   
136   Int_t nSkippedSeeds = 0;
137   
138   // loop over esd tracks
139   for (Int_t t=0; t<nTracks; t++)
140   {
141
142     AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (fESD->GetTrack(t));
143     if (!esdTrack)
144     {
145       AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", t));
146       continue;
147     }
148     
149     AliESDfriendTrack* friendtrack = const_cast<AliESDfriendTrack*> (dynamic_cast<const AliESDfriendTrack*> (esdTrack->GetFriendTrack()));
150     if (!friendtrack)
151     {
152       AliDebug(AliLog::kError, Form("ERROR: Could not retrieve friend of track %d.", t));
153       continue;
154     }
155     
156     const AliTPCseed* seed = dynamic_cast<const AliTPCseed*> (friendtrack->GetCalibObject(0));
157     if (!seed)
158     {
159       AliDebug(AliLog::kDebug, Form("ERROR: Could not retrieve seed of track %d.", t));
160       nSkippedSeeds++;
161       continue;
162     }
163     
164     for (Int_t clusterID = 0; clusterID < 160; clusterID++)
165     {
166       AliTPCclusterMI* cluster = seed->GetClusterPointer(clusterID);
167       if (!cluster)
168       {
169         //AliDebug(AliLog::kError, Form("ERROR: Could not retrieve cluster %d of track %d.", clusterID, t));
170         continue;
171       }
172       
173       //AliDebug(AliLog::kDebug, Form("We found a cluster from sector %d", cluster->GetDetector()));
174
175       Int_t detector = cluster->GetDetector();
176       
177       if (detector < 0 || detector >= kTPCSectors) {
178         AliDebug(AliLog::kDebug, Form("We found a cluster from invalid sector %d", detector));
179         continue;
180       }
181       
182       // TODO: find a clever way to handle the time      
183       //      if (fESD->GetTimeStamp()<1160000000)
184       //        continue;
185
186       if (!fClusterHistograms[detector])
187       {
188         fClusterHistograms[detector] = new AliTPCClusterHistograms(detector,"",fESD->GetTimeStamp(),fESD->GetTimeStamp()+7*60*60);
189       }
190       
191       fClusterHistograms[detector]->FillCluster(cluster, fESD->GetTimeStamp());
192     }
193   }
194   
195   if (nSkippedSeeds > 0)
196     printf("WARNING: The seed was not found for %d out of %d tracks.\n", nSkippedSeeds, nTracks);
197    
198   return kTRUE;
199 }
200
201 void AliROCESDAnalysisSelector::SlaveTerminate()
202 {
203   //
204   
205   for (Int_t i=0; i<kTPCSectors; i++)
206     if (fClusterHistograms[i])
207         fOutput->Add(fClusterHistograms[i]);
208
209
210
211 void AliROCESDAnalysisSelector::Terminate()
212 {
213   // TODO read from output list for PROOF
214     
215   TNamed* comment = dynamic_cast<TNamed*>(fTree->GetUserInfo()->FindObject("comment"));
216
217   if (comment)
218     AliDebug(AliLog::kInfo, Form("INFO: Found comment in input list: %s \n", comment->GetTitle()));
219
220   TFile* file = TFile::Open(Form("rocESD_%s.root",comment->GetTitle()), "RECREATE");
221   
222   for (Int_t i=0; i<kTPCSectors; i++)
223     if (fClusterHistograms[i]) {
224       fClusterHistograms[i]->SaveHistograms();
225       TCanvas* c = fClusterHistograms[i]->DrawHistograms();
226       c->SaveAs(Form("plots_%s_%s.eps",comment->GetTitle(),c->GetName()));
227       c->SaveAs(Form("plots_%s_%s.gif",comment->GetTitle(),c->GetName()));
228
229       c->Close();
230       delete c;
231     }
232   file->Close();
233