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