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.
36 class AliTrackReference;
42 class AliFMDRawReader;
58 //___________________________________________________________________
59 /** @class AliFMDInput
60 @brief Base class for reading in various FMD data.
61 The class loops over all found events. For each event the
62 specified data is read in. The class then loops over all
63 elements of the read data, and process these with user defined
66 struct DigitInput : public AliFMDInput
73 fHist = new TH1F("adc", "ADC spectra", 1024, -.5, 1023.5);
76 Bool_t ProcessDigit(AliFMDDigit* d)
78 fHist->Fill(d->Counts());
81 // After processing all events, display spectrum
95 This class allows for writing small scripts, that can be compiled
96 with AcLIC, to do all sorts of tests, quick prototyping, and so
97 on. It has proven to be quiet useful. One can load more than
98 one type of data in one derived class, to for example to make
99 comparisons between hits and reconstructed points. See also the
100 various scripts in @c FMD/scripts.
103 class AliFMDInput : public TNamed
106 /** The kinds of data that can be read in. */
109 kKinematics, // Kinematics (from sim)
111 kSDigits, // Summable digits
112 kHeader, // Header information
113 kRecPoints, // Reconstructed points
115 kRaw, // Read raw data
116 kGeometry, // Not really a tree
117 kTracks, // Hits and tracs - for BG study
118 kTrackRefs, // Track references - also for BG study
119 kRawCalib, // Read raws and calibrate them
125 @param gAliceFile galice file */
126 AliFMDInput(const char* gAliceFile);
128 virtual ~AliFMDInput() {}
130 /** Add a data type to load
131 @param tree Data to load */
132 virtual void AddLoad(ETrees tree) { SETBIT(fTreeMask, tree); }
133 /** Remove a data type to load
134 @param tree Data to @e not load */
135 virtual void RemoveLoad(ETrees tree) { CLRBIT(fTreeMask, tree); }
136 /** @return # of available events */
137 virtual Int_t NEvents() const;
139 /** Initialize the class. If a user class overloads this member
140 function, then this @e must be explicitly called
141 @return @c false on error */
142 virtual Bool_t Init();
143 /** Callled at the beginning of each event. If a user class
144 overloads this member function, then this @e must be explicitly
146 @param event Event number
147 @return @c false on error */
148 virtual Bool_t Begin(Int_t event);
149 /** Process one event. This loops over all the loaded data. Users
150 can overload this member function, but then it's @e strongly
151 recommended to explicitly call this classes version.
152 @return @c false on error */
153 virtual Bool_t Event();
154 /** Called at the end of each event.
155 @return @c false on error */
156 virtual Bool_t End();
157 /** Called at the end of the run.
158 @return @c false on error */
159 virtual Bool_t Finish() { return kTRUE; }
161 @return @c false on error */
162 virtual Bool_t Run();
164 /** Loop over all hits, and call ProcessHit with that hit, and
165 optionally the corresponding kinematics track.
166 @return @c false on error */
167 virtual Bool_t ProcessHits();
168 /** Loop over all track refs, and call ProcessTrackRef with that hit, and
169 optionally the corresponding kinematics track.
170 @return @c false on error */
171 virtual Bool_t ProcessTrackRefs();
172 /** Loop over all tracks, and call ProcessTrack with each hit for
174 @return @c false on error */
175 virtual Bool_t ProcessTracks();
176 /** Loop over all digits, and call ProcessDigit for each digit.
177 @return @c false on error */
178 virtual Bool_t ProcessDigits();
179 /** Loop over all summable digits, and call ProcessSDigit for each
181 @return @c false on error */
182 virtual Bool_t ProcessSDigits();
183 /** Loop over all digits read from raw data files, and call
184 ProcessRawDigit for each digit.
185 @return @c false on error */
186 virtual Bool_t ProcessRawDigits();
187 /** Loop over all digits read from raw data files, and call
188 ProcessRawDigit for each digit.
189 @return @c false on error */
190 virtual Bool_t ProcessRawCalibDigits();
191 /** Loop over all reconstructed points, and call ProcessRecPoint for
192 each reconstructed point.
193 @return @c false on error */
194 virtual Bool_t ProcessRecPoints();
195 /** Loop over all ESD data, and call ProcessESD for each entry.
196 @return @c false on error */
197 virtual Bool_t ProcessESDs();
198 /** Loop over all strips and ask user routine to supply the data.
199 @return @c false on error */
200 virtual Bool_t ProcessUsers();
202 /** Process one hit, and optionally it's corresponding kinematics
203 track. Users should over this to process each hit.
205 @param p Associated track
206 @return @c false on error */
207 virtual Bool_t ProcessHit(AliFMDHit* h, TParticle* p);
208 /** Process one track reference, and optionally it's corresponding kinematics
209 track. Users should overload this to process each track reference.
210 @param trackRef Track Reference
211 @param track Associated track
212 @return @c false on error */
213 virtual Bool_t ProcessTrackRef(AliTrackReference* trackRef, TParticle* track);
214 /** Process one hit per track. Users should over this to process
216 @param i Track number
218 @param h Associated Hit
219 @return @c false on error */
220 virtual Bool_t ProcessTrack(Int_t i, TParticle* p, AliFMDHit* h);
221 /** Process one digit. Users should over this to process each
224 @return @c false on error */
225 virtual Bool_t ProcessDigit(AliFMDDigit* digit);
226 /** Process one summable digit. Users should over this to process
228 @param sdigit Summable digit
229 @return @c false on error */
230 virtual Bool_t ProcessSDigit(AliFMDSDigit* sdigit);
231 /** Process one digit from raw data files. Users should over this
232 to process each raw digit.
233 @param digit Raw digit
234 @return @c false on error */
235 virtual Bool_t ProcessRawDigit(AliFMDDigit* digit);
236 /** Process one digit from raw data files. Users should over this
237 to process each raw digit.
238 @param digit Raw digit
239 @return @c false on error */
240 virtual Bool_t ProcessRawCalibDigit(AliFMDDigit* digit);
241 /** Process one reconstructed point. Users should over this to
242 process each reconstructed point.
243 @param point Reconstructed point
244 @return @c false on error */
245 virtual Bool_t ProcessRecPoint(AliFMDRecPoint* point);
246 /** Process ESD data for the FMD. Users should overload this to
248 @param d Detector number (1-3)
249 @param r Ring identifier ('I' or 'O')
250 @param s Sector number (0-19, or 0-39)
251 @param t Strip number (0-511, or 0-255)
252 @param eta Psuedo-rapidity
253 @param mult Psuedo-multiplicity
254 @return @c false on error */
255 virtual Bool_t ProcessESD(UShort_t d, Char_t r, UShort_t s, UShort_t t,
256 Float_t eta, Float_t mult);
257 /** Process User data for the FMD. Users should overload this to
259 @param d Detector number (1-3)
260 @param r Ring identifier ('I' or 'O')
261 @param s Sector number (0-19, or 0-39)
262 @param t Strip number (0-511, or 0-255)
264 @return @c false on error */
265 virtual Bool_t ProcessUser(UShort_t d, Char_t r, UShort_t s, UShort_t t,
267 /** Service function to make a logarithmic axis.
268 @param n Number of bins
269 @param min Minimum of axis
270 @param max Maximum of axis.
271 @return An array with the bin boundaries. */
272 static TArrayF MakeLogScale(Int_t n, Double_t min, Double_t max);
274 /** Set the raw data input
275 @param file File name - if empty, assume simulated raw. */
276 void SetRawFile(const char* file) { if (file) fRawFile = file; }
280 @param o Object to copy from */
281 AliFMDInput(const AliFMDInput& o)
315 /** Assignement operator
316 @return REference to this */
317 AliFMDInput& operator=(const AliFMDInput&) { return *this; }
319 * Get user supplued data
328 virtual Float_t GetSignal(UShort_t d, Char_t r, UShort_t s, UShort_t t);
330 TString fGAliceFile; // File name of gAlice file
331 AliRunLoader* fLoader; // Loader of FMD data
332 AliRun* fRun; // Run information
333 AliStack* fStack; // Stack of particles
334 AliLoader* fFMDLoader; // Loader of FMD data
335 AliRawReader* fReader; // Raw data reader
336 AliFMDRawReader* fFMDReader; // FMD raw reader
337 AliFMD* fFMD; // FMD object
338 AliESDFMD* fESD; // FMD ESD data
339 AliESDEvent* fESDEvent; // ESD Event object.
340 TTree* fTreeE; // Header tree
341 TTree* fTreeH; // Hits tree
342 TTree* fTreeTR; // Track Reference tree
343 TTree* fTreeD; // Digit tree
344 TTree* fTreeS; // SDigit tree
345 TTree* fTreeR; // RecPoint tree
346 TTree* fTreeA; // Raw data tree
347 TChain* fChainE; // Chain of ESD's
348 TClonesArray* fArrayE; // Event info array
349 TClonesArray* fArrayH; // Hit info array
350 TClonesArray* fArrayTR; // Hit info array
351 TClonesArray* fArrayD; // Digit info array
352 TClonesArray* fArrayS; // SDigit info array
353 TClonesArray* fArrayR; // Rec points info array
354 TClonesArray* fArrayA; // Raw data (digits) info array
355 AliHeader* fHeader; // Header
356 TGeoManager* fGeoManager; // Geometry manager
357 Int_t fTreeMask; // Which tree's to load
358 TString fRawFile; // Raw input file
359 Bool_t fIsInit; // Have we been initialized
360 Int_t fEventCount; // Event counter
361 ClassDef(AliFMDInput,0) //Hits for detector FMD
364 inline Bool_t AliFMDInput::ProcessHit(AliFMDHit*,TParticle*) { return kTRUE; }
365 inline Bool_t AliFMDInput::ProcessTrackRef(AliTrackReference*,
366 TParticle*) { return kTRUE; }
367 inline Bool_t AliFMDInput::ProcessTrack(Int_t,TParticle*,
368 AliFMDHit*) { return kTRUE; }
369 inline Bool_t AliFMDInput::ProcessDigit(AliFMDDigit*) { return kTRUE; }
370 inline Bool_t AliFMDInput::ProcessSDigit(AliFMDSDigit*) { return kTRUE; }
371 inline Bool_t AliFMDInput::ProcessRawDigit(AliFMDDigit*) { return kTRUE; }
372 inline Bool_t AliFMDInput::ProcessRawCalibDigit(AliFMDDigit*) { return kTRUE; }
373 inline Bool_t AliFMDInput::ProcessRecPoint(AliFMDRecPoint*) { return kTRUE; }
374 inline Bool_t AliFMDInput::ProcessESD(UShort_t,Char_t,UShort_t,UShort_t,
375 Float_t,Float_t) { return kTRUE; }
376 inline Bool_t AliFMDInput::ProcessUser(UShort_t,Char_t,UShort_t,UShort_t,
377 Float_t) { return kTRUE; }
378 inline Float_t AliFMDInput::GetSignal(UShort_t, Char_t, UShort_t, UShort_t) {
383 //____________________________________________________________________