Added documentation of each file.
[u/mrichter/AliRoot.git] / FMD / AliFMDInput.h
1 #ifndef AliFMDInput_H
2 #define AliFMDInput_H
3 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights
4  * reserved. 
5  *
6  * See cxx source for full Copyright notice                               
7  */
8 /** @file    AliFMDInput.h
9     @author  Christian Holm Christensen <cholm@nbi.dk>
10     @date    Mon Mar 27 12:42:40 2006
11     @brief   FMD utility classes for reading FMD data
12 */
13 //___________________________________________________________________
14 /** @defgroup FMD_util Utility classes. 
15
16     The classes defined here, are utility classes for reading in data
17     for the FMD.  They are put in a seperate library to not polute the
18     normal libraries.  The classes are intended to be used as base
19     classes for customized class that do some sort of analysis on the
20     various types of data produced by the FMD.
21 */
22 #include <TObject.h>
23 #ifndef ROOT_TString
24 # include <TString.h>
25 #endif
26 class AliRunLoader;
27 class AliLoader;
28 class AliStack;
29 class AliRun;
30 class AliRawReader;
31 class AliFMD;
32 class AliFMDHit;
33 class AliFMDDigit;
34 class AliFMDSDigit;
35 class AliFMDRecPoint;
36 class AliESD;
37 class AliESDFMD;
38 class TString;
39 class TClonesArray;
40 class TTree;
41 class TGeoManager;
42 class TParticle;
43 class TChain;
44
45 //___________________________________________________________________
46 /** @class AliFMDInput 
47     @brief Base class for reading in various FMD data. 
48     The class loops over all found events.  For each event the
49     specified data is read in.   The class then loops over all
50     elements of the read data, and process these with user defined
51     code. 
52     @code 
53     struct DigitInput : public AliFMDInput 
54     { 
55       DigitInput() 
56       { 
57          // Load digits 
58          AddLoad(kDigits);
59          // Make a histogram 
60          fHist = new TH1F("adc", "ADC spectra", 1024, -.5, 1023.5);
61       }
62       // Process one digit. 
63       Bool_t ProcessDigit(AliFMDDigit* d)
64       { 
65          fHist->Fill(d->Counts());
66          return kTRUE;
67       }
68       // After processing all events, display spectrum
69       Bool_t Finish()
70       {
71         fHist->Draw();
72       }
73       TH1F* fHist;
74     };
75
76     void AdcSpectrum()
77     { 
78        DigitInput di;
79        di.Run();
80     }
81     @endcode
82     This class allows for writing small scripts, that can be compiled
83     with AcLIC, to do all sorts of tests, quick prototyping, and so
84     on.  It has proven to be quiet useful.   One can load more than
85     one type of data in one derived class, to for example to make
86     comparisons between hits and reconstructed points.   See also the
87     various scripts in @c FMD/scripts. 
88     @ingroup FMD_util
89  */
90 class AliFMDInput : public TObject
91 {
92 public:
93   /** The kinds of data that can be read in. */
94   enum ETrees {
95     kHits       = 1,  // Hits
96     kKinematics,      // Kinematics (from sim)
97     kDigits,          // Digits
98     kSDigits,         // Summable digits 
99     kHeader,          // Header information 
100     kRecPoints,       // Reconstructed points
101     kESD,             // Load ESD's
102     kRaw,             // Read raw data 
103     kGeometry         // Not really a tree 
104   };
105   /** CTOR  */
106   AliFMDInput();
107   /** CTOR
108       @param gAliceFile galice file  */
109   AliFMDInput(const char* gAliceFile);
110   /** DTOR */
111   virtual ~AliFMDInput() {}
112
113   /** Add a data type to load 
114       @param tree Data to load */
115   virtual void   AddLoad(ETrees tree)     { SETBIT(fTreeMask, tree); }
116   /** Remove a data type to load 
117       @param tree Data to @e not load */
118   virtual void   RemoveLoad(ETrees tree)  { CLRBIT(fTreeMask, tree); }
119   /** @return # of available events */
120   virtual Int_t  NEvents() const;
121
122   /** Initialize the class.  If a user class overloads this member
123       function, then this @e must be explicitly called
124       @return @c false on error */
125   virtual Bool_t Init();
126   /** Callled at the beginning of each event. If a user class
127       overloads this member  function, then this @e must be explicitly
128       called. 
129       @param event Event number
130       @return @c false on error */
131   virtual Bool_t Begin(Int_t event);
132   /** Process one event.  This loops over all the loaded data.   Users
133       can overload this member function, but then it's @e strongly
134       recommended to explicitly call this classes version. 
135       @return @c false on error  */
136   virtual Bool_t Event();
137   /** Called at the end of each event. 
138       @return @c false on error  */
139   virtual Bool_t End();
140   /** Called at the end of the run.
141       @return  @c false on error */
142   virtual Bool_t Finish() { return kTRUE; }
143   /** Run a full job. 
144       @return @c false on error  */
145   virtual Bool_t Run();
146
147   /** Loop over all hits, and call ProcessHit with that hit, and
148       optionally the corresponding kinematics track. 
149       @return @c false on error  */
150   virtual Bool_t ProcessHits();
151   /** Loop over all digits, and call ProcessDigit for each digit.
152       @return @c false on error  */
153   virtual Bool_t ProcessDigits();
154   /** Loop over all summable digits, and call ProcessSDigit for each
155       digit. 
156       @return @c false on error  */
157   virtual Bool_t ProcessSDigits();
158   /** Loop over all digits read from raw data files, and call
159       ProcessRawDigit for each digit. 
160       @return @c false on error  */
161   virtual Bool_t ProcessRawDigits();
162   /** Loop over all reconstructed points, and call ProcessRecPoint for
163       each reconstructed point. 
164       @return @c false on error  */
165   virtual Bool_t ProcessRecPoints();
166
167   /** Process one hit, and optionally it's corresponding kinematics
168       track.  Users should over this to process each hit. 
169       @return  @c false on error   */
170   virtual Bool_t ProcessHit(AliFMDHit*, TParticle*)  { return kTRUE; }
171   /** Process one digit.  Users should over this to process each digit. 
172       @return  @c false on error   */
173   virtual Bool_t ProcessDigit(AliFMDDigit*)          { return kTRUE; }
174   /** Process one summable digit.  Users should over this to process
175       each summable digit.  
176       @return  @c false on error   */
177   virtual Bool_t ProcessSDigit(AliFMDSDigit*)        { return kTRUE; }
178   /** Process one digit from raw data files.  Users should over this
179       to process each raw digit.  
180       @return  @c false on error   */
181   virtual Bool_t ProcessRawDigit(AliFMDDigit*)       { return kTRUE; }
182   /** Process one reconstructed point.  Users should over this to
183       process each reconstructed point.  
184       @return  @c false on error   */
185   virtual Bool_t ProcessRecPoint(AliFMDRecPoint*)    { return kTRUE; }
186   /** Process ESD data for the FMD.  Users should overload this to
187       deal with ESD data. 
188       @return  @c false on error  */
189   virtual Bool_t ProcessESD(AliESDFMD*)              { return kTRUE; }
190   
191 protected:
192   TString       fGAliceFile; // File name of gAlice file
193   AliRunLoader* fLoader;     // Loader of FMD data 
194   AliRun*       fRun;        // Run information
195   AliStack*     fStack;      // Stack of particles 
196   AliLoader*    fFMDLoader;  // Loader of FMD data 
197   AliRawReader* fReader;     // Raw data reader 
198   AliFMD*       fFMD;        // FMD object
199   AliESD*       fMainESD;    // ESD Object
200   AliESDFMD*    fESD;        // FMD ESD data  
201   TTree*        fTreeE;      // Header tree 
202   TTree*        fTreeH;      // Hits tree
203   TTree*        fTreeD;      // Digit tree 
204   TTree*        fTreeS;      // SDigit tree 
205   TTree*        fTreeR;      // RecPoint tree
206   TTree*        fTreeA;      // Raw data tree
207   TChain*       fChainE;     // Chain of ESD's
208   TClonesArray* fArrayE;     // Event info array
209   TClonesArray* fArrayH;     // Hit info array
210   TClonesArray* fArrayD;     // Digit info array
211   TClonesArray* fArrayS;     // SDigit info array
212   TClonesArray* fArrayR;     // Rec points info array
213   TClonesArray* fArrayA;     // Raw data (digits) info array
214   TGeoManager*  fGeoManager; // Geometry manager
215   Int_t         fTreeMask;   // Which tree's to load
216   Bool_t        fIsInit;
217   ClassDef(AliFMDInput,0)  //Hits for detector FMD
218 };
219
220
221 //____________________________________________________________________
222 class AliFMDHit;
223 /** @brief Class to read FMD hits
224  */
225 class AliFMDInputHits : public AliFMDInput 
226 {
227 public:
228   /** Constructor
229       @param file Name of @c gAlice file */
230   AliFMDInputHits(const char* file="galice.root") 
231     : AliFMDInput(file) { AddLoad(kHits); }
232   ClassDef(AliFMDInputHits, 0);
233 };
234
235 //____________________________________________________________________
236 class AliFMDDigit;
237 /** @brief Class to read FMD digits
238  */
239 class AliFMDInputDigits : public AliFMDInput 
240 {
241 public:
242   /** Constructor
243       @param file Name of @c gAlice file */
244   AliFMDInputDigits(const char* file="galice.root")
245     : AliFMDInput(file) { AddLoad(kDigits); }
246   ClassDef(AliFMDInputDigits, 0);
247 };
248
249 //____________________________________________________________________
250 class AliFMDSDigit;
251 /** @brief Class to read FMD summable digits
252  */
253 class AliFMDInputSDigits : public AliFMDInput 
254 {
255 public:
256   /** Constructor
257       @param file Name of @c gAlice file */
258   AliFMDInputSDigits(const char* file="galice.root") 
259     : AliFMDInput(file) { AddLoad(kSDigits); }
260   ClassDef(AliFMDInputSDigits, 0);
261 };
262
263 //____________________________________________________________________
264 /** @brief Class to read FMD raw data
265  */
266 class AliFMDInputRaw : public AliFMDInput 
267 {
268 public:
269   /** Constructor
270       @param file Name of @c gAlice file */
271   AliFMDInputRaw(const char* file="galice.root") 
272     : AliFMDInput(file) { AddLoad(kRaw); }
273   ClassDef(AliFMDInputRaw, 0);
274 };
275
276 //____________________________________________________________________
277 class AliFMDRecPoint;
278 /** @brief Class to read FMD reconstructed data
279  */
280 class AliFMDInputRecPoints : public AliFMDInput 
281 {
282 public:
283   /** Constructor
284       @param file Name of @c gAlice file */
285   AliFMDInputRecPoints(const char* file="galice.root") 
286     : AliFMDInput(file) { AddLoad(kRecPoints); }
287   ClassDef(AliFMDInputRecPoints, 0);
288 };
289
290 #endif
291 //____________________________________________________________________
292 //
293 // Local Variables:
294 //   mode: C++
295 // End:
296 //
297 // EOF
298 //