3 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights
6 * See cxx source for full Copyright notice
8 //___________________________________________________________________
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.
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
20 //___________________________________________________________________
21 /** @defgroup FMD_util Utility classes.
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.
55 //___________________________________________________________________
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
63 struct DigitInput : public AliFMDInput
70 fHist = new TH1F("adc", "ADC spectra", 1024, -.5, 1023.5);
73 Bool_t ProcessDigit(AliFMDDigit* d)
75 fHist->Fill(d->Counts());
78 // After processing all events, display spectrum
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.
100 class AliFMDInput : public TObject
103 /** The kinds of data that can be read in. */
106 kKinematics, // Kinematics (from sim)
108 kSDigits, // Summable digits
109 kHeader, // Header information
110 kRecPoints, // Reconstructed points
112 kRaw, // Read raw data
113 kGeometry, // Not really a tree
114 kTracks // Hits and tracs - for BG study
119 @param gAliceFile galice file */
120 AliFMDInput(const char* gAliceFile);
122 virtual ~AliFMDInput() {}
124 /** Add a data type to load
125 @param tree Data to load */
126 virtual void AddLoad(ETrees tree) { SETBIT(fTreeMask, tree); }
127 /** Remove a data type to load
128 @param tree Data to @e not load */
129 virtual void RemoveLoad(ETrees tree) { CLRBIT(fTreeMask, tree); }
130 /** @return # of available events */
131 virtual Int_t NEvents() const;
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 */
136 virtual Bool_t Init();
137 /** Callled at the beginning of each event. If a user class
138 overloads this member function, then this @e must be explicitly
140 @param event Event number
141 @return @c false on error */
142 virtual Bool_t Begin(Int_t event);
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 */
147 virtual Bool_t Event();
148 /** Called at the end of each event.
149 @return @c false on error */
150 virtual Bool_t End();
151 /** Called at the end of the run.
152 @return @c false on error */
153 virtual Bool_t Finish() { return kTRUE; }
155 @return @c false on error */
156 virtual Bool_t Run();
158 /** Loop over all hits, and call ProcessHit with that hit, and
159 optionally the corresponding kinematics track.
160 @return @c false on error */
161 virtual Bool_t ProcessHits();
162 /** Loop over all tracks, and call ProcessTrack with each hit for
164 @return @c false on error */
165 virtual Bool_t ProcessTracks();
166 /** Loop over all digits, and call ProcessDigit for each digit.
167 @return @c false on error */
168 virtual Bool_t ProcessDigits();
169 /** Loop over all summable digits, and call ProcessSDigit for each
171 @return @c false on error */
172 virtual Bool_t ProcessSDigits();
173 /** Loop over all digits read from raw data files, and call
174 ProcessRawDigit for each digit.
175 @return @c false on error */
176 virtual Bool_t ProcessRawDigits();
177 /** Loop over all reconstructed points, and call ProcessRecPoint for
178 each reconstructed point.
179 @return @c false on error */
180 virtual Bool_t ProcessRecPoints();
181 /** Loop over all ESD data, and call ProcessESD for each entry.
182 @return @c false on error */
183 virtual Bool_t ProcessESDs();
185 /** Process one hit, and optionally it's corresponding kinematics
186 track. Users should over this to process each hit.
188 @param p Associated track
189 @return @c false on error */
190 virtual Bool_t ProcessHit(AliFMDHit* h, TParticle* p);
191 /** Process one hit per track. Users should over this to process
193 @param i Track number
195 @param h Associated Hit
196 @return @c false on error */
197 virtual Bool_t ProcessTrack(Int_t i, TParticle* p, AliFMDHit* h);
198 /** Process one digit. Users should over this to process each
201 @return @c false on error */
202 virtual Bool_t ProcessDigit(AliFMDDigit* digit);
203 /** Process one summable digit. Users should over this to process
205 @param sdigit Summable digit
206 @return @c false on error */
207 virtual Bool_t ProcessSDigit(AliFMDSDigit* sdigit);
208 /** Process one digit from raw data files. Users should over this
209 to process each raw digit.
210 @param digit Raw digit
211 @return @c false on error */
212 virtual Bool_t ProcessRawDigit(AliFMDDigit* digit);
213 /** Process one reconstructed point. Users should over this to
214 process each reconstructed point.
215 @param point Reconstructed point
216 @return @c false on error */
217 virtual Bool_t ProcessRecPoint(AliFMDRecPoint* point);
218 /** Process ESD data for the FMD. Users should overload this to
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
226 @return @c false on error */
227 virtual Bool_t ProcessESD(UShort_t, Char_t, UShort_t, UShort_t,
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);
237 @param o Object to copy from */
238 AliFMDInput(const AliFMDInput& o)
266 /** Assignement operator
267 @return REference to this */
268 AliFMDInput& operator=(const AliFMDInput&) { return *this; }
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
275 AliRawReader* fReader; // Raw data reader
276 AliFMD* fFMD; // FMD object
277 AliESDFMD* fESD; // FMD ESD data
278 AliESDEvent* fESDEvent; // ESD Event object.
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
284 TTree* fTreeA; // Raw data tree
285 TChain* fChainE; // Chain of ESD's
286 TClonesArray* fArrayE; // Event info array
287 TClonesArray* fArrayH; // Hit info array
288 TClonesArray* fArrayD; // Digit info array
289 TClonesArray* fArrayS; // SDigit info array
290 TClonesArray* fArrayR; // Rec points info array
291 TClonesArray* fArrayA; // Raw data (digits) info array
292 TGeoManager* fGeoManager; // Geometry manager
293 Int_t fTreeMask; // Which tree's to load
294 Bool_t fIsInit; // Have we been initialized
295 Int_t fEventCount; // Event counter
296 ClassDef(AliFMDInput,0) //Hits for detector FMD
299 inline Bool_t AliFMDInput::ProcessHit(AliFMDHit*,TParticle*) { return kTRUE; }
300 inline Bool_t AliFMDInput::ProcessTrack(Int_t,TParticle*,
301 AliFMDHit*) { return kTRUE; }
302 inline Bool_t AliFMDInput::ProcessDigit(AliFMDDigit*) { return kTRUE; }
303 inline Bool_t AliFMDInput::ProcessSDigit(AliFMDSDigit*) { return kTRUE; }
304 inline Bool_t AliFMDInput::ProcessRawDigit(AliFMDDigit*) { return kTRUE; }
305 inline Bool_t AliFMDInput::ProcessRecPoint(AliFMDRecPoint*) { return kTRUE; }
306 inline Bool_t AliFMDInput::ProcessESD(UShort_t,Char_t,UShort_t,UShort_t,
307 Float_t,Float_t) { return kTRUE; }
311 //____________________________________________________________________