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;
57 class TSystemDirectory;
59 //___________________________________________________________________
60 /** @class AliFMDInput
61 @brief Base class for reading in various FMD data.
62 The class loops over all found events. For each event the
63 specified data is read in. The class then loops over all
64 elements of the read data, and process these with user defined
67 struct DigitInput : public AliFMDInput
74 fHist = new TH1F("adc", "ADC spectra", 1024, -.5, 1023.5);
77 Bool_t ProcessDigit(AliFMDDigit* d)
79 fHist->Fill(d->Counts());
82 // After processing all events, display spectrum
96 This class allows for writing small scripts, that can be compiled
97 with AcLIC, to do all sorts of tests, quick prototyping, and so
98 on. It has proven to be quiet useful. One can load more than
99 one type of data in one derived class, to for example to make
100 comparisons between hits and reconstructed points. See also the
101 various scripts in @c FMD/scripts.
104 class AliFMDInput : public TNamed
107 /** The kinds of data that can be read in. */
110 kKinematics, // Kinematics (from sim)
112 kSDigits, // Summable digits
113 kHeader, // Header information
114 kRecPoints, // Reconstructed points
116 kRaw, // Read raw data
117 kGeometry, // Not really a tree
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;
138 /** @return true if passed tree is loaded */
139 virtual Bool_t IsLoaded(ETrees tree)const { return TESTBIT(fTreeMask, tree); }
141 * Set the trees to load.
143 * @param mask Bit mask of trees to load. Should be constructed
147 * UInt_t mask = ((1 << AliFMDInput::kHits) |
148 * (1 << AliFMDInput::kDigits) |
149 * (1 << AliFMDInput::kSDigits));
152 virtual void SetLoads(UInt_t mask);
154 * Set the trees to load.
156 * @param mask A comma or space separated list of trees to load.
157 * The case is not important, and a short from of the tree name can
160 virtual void SetLoads(const char* mask);
162 * Get a string that represents the loaded trees
164 * @param dataOnly If true, then only show data
166 * @return String representing the loaded trees.
168 virtual const char* LoadedString(Bool_t dataOnly=false) const;
170 /** Initialize the class. If a user class overloads this member
171 function, then this @e must be explicitly called
172 @return @c false on error */
173 virtual Bool_t Init();
174 /** Callled at the beginning of each event. If a user class
175 overloads this member function, then this @e must be explicitly
177 @param event Event number
178 @return @c false on error */
179 virtual Bool_t Begin(Int_t event);
180 /** Process one event. This loops over all the loaded data. Users
181 can overload this member function, but then it's @e strongly
182 recommended to explicitly call this classes version.
183 @return @c false on error */
184 virtual Bool_t Event();
185 /** Called at the end of each event.
186 @return @c false on error */
187 virtual Bool_t End();
188 /** Called at the end of the run.
189 @return @c false on error */
190 virtual Bool_t Finish() { return kTRUE; }
192 @return @c false on error */
193 virtual Bool_t Run(UInt_t maxEvents=0);
195 /** Loop over all hits, and call ProcessHit with that hit, and
196 optionally the corresponding kinematics track.
197 @return @c false on error */
198 virtual Bool_t ProcessHits();
199 /** Loop over all track refs, and call ProcessTrackRef with that hit, and
200 optionally the corresponding kinematics track.
201 @return @c false on error */
202 virtual Bool_t ProcessTrackRefs();
203 /** Loop over all tracks, and call ProcessTrack with each hit for
205 @return @c false on error */
206 virtual Bool_t ProcessTracks();
207 /** Loop over all tracks, and call ProcessTrack with each hit for
209 @return @c false on error */
210 virtual Bool_t ProcessStack();
211 /** Loop over all digits, and call ProcessDigit for each digit.
212 @return @c false on error */
213 virtual Bool_t ProcessDigits();
214 /** Loop over all summable digits, and call ProcessSDigit for each
216 @return @c false on error */
217 virtual Bool_t ProcessSDigits();
218 /** Loop over all digits read from raw data files, and call
219 ProcessRawDigit for each digit.
220 @return @c false on error */
221 virtual Bool_t ProcessRawDigits();
222 /** Loop over all digits read from raw data files, and call
223 ProcessRawDigit for each digit.
224 @return @c false on error */
225 virtual Bool_t ProcessRawCalibDigits();
226 /** Loop over all reconstructed points, and call ProcessRecPoint for
227 each reconstructed point.
228 @return @c false on error */
229 virtual Bool_t ProcessRecPoints();
230 /** Loop over all ESD data, and call ProcessESD for each entry.
231 @return @c false on error */
232 virtual Bool_t ProcessESDs();
233 /** Loop over all strips and ask user routine to supply the data.
234 @return @c false on error */
235 virtual Bool_t ProcessUsers();
237 /** Process one hit, and optionally it's corresponding kinematics
238 track. Users should over this to process each hit.
240 @param p Associated track
241 @return @c false on error */
242 virtual Bool_t ProcessHit(AliFMDHit* h, TParticle* p);
243 /** Process one track reference, and optionally it's corresponding kinematics
244 track. Users should overload this to process each track reference.
245 @param trackRef Track Reference
246 @param track Associated track
247 @return @c false on error */
248 virtual Bool_t ProcessTrackRef(AliTrackReference* trackRef, TParticle* track);
249 /** Process one hit per track. Users should over this to process
251 @param i Track number
253 @param h Associated Hit
254 @return @c false on error */
255 virtual Bool_t ProcessTrack(Int_t i, TParticle* p, AliFMDHit* h);
256 /** Process stack particle
257 @param i Track number
259 @return @c false on error */
260 virtual Bool_t ProcessParticle(Int_t i , TParticle* p);
261 /** Process one digit. Users should over this to process each
264 @return @c false on error */
265 virtual Bool_t ProcessDigit(AliFMDDigit* digit);
266 /** Process one summable digit. Users should over this to process
268 @param sdigit Summable digit
269 @return @c false on error */
270 virtual Bool_t ProcessSDigit(AliFMDSDigit* sdigit);
271 /** Process one digit from raw data files. Users should over this
272 to process each raw digit.
273 @param digit Raw digit
274 @return @c false on error */
275 virtual Bool_t ProcessRawDigit(AliFMDDigit* digit);
276 /** Process one digit from raw data files. Users should over this
277 to process each raw digit.
278 @param digit Raw digit
279 @return @c false on error */
280 virtual Bool_t ProcessRawCalibDigit(AliFMDDigit* digit);
281 /** Process one reconstructed point. Users should over this to
282 process each reconstructed point.
283 @param point Reconstructed point
284 @return @c false on error */
285 virtual Bool_t ProcessRecPoint(AliFMDRecPoint* point);
286 /** Process ESD data for the FMD. Users should overload this to
288 @param d Detector number (1-3)
289 @param r Ring identifier ('I' or 'O')
290 @param s Sector number (0-19, or 0-39)
291 @param t Strip number (0-511, or 0-255)
292 @param eta Psuedo-rapidity
293 @param mult Psuedo-multiplicity
294 @return @c false on error */
295 virtual Bool_t ProcessESD(UShort_t d, Char_t r, UShort_t s, UShort_t t,
296 Float_t eta, Float_t mult);
297 /** Process User data for the FMD. Users should overload this to
299 @param d Detector number (1-3)
300 @param r Ring identifier ('I' or 'O')
301 @param s Sector number (0-19, or 0-39)
302 @param t Strip number (0-511, or 0-255)
304 @return @c false on error */
305 virtual Bool_t ProcessUser(UShort_t d, Char_t r, UShort_t s, UShort_t t,
307 /** Service function to make a logarithmic axis.
308 @param n Number of bins
309 @param min Minimum of axis
310 @param max Maximum of axis.
311 @return An array with the bin boundaries. */
312 static TArrayF MakeLogScale(Int_t n, Double_t min, Double_t max);
314 /** Set the raw data input
315 @param file File name - if empty, assume simulated raw. */
316 void SetRawFile(const char* file) { if (file) fRawFile = file; }
317 void SetInputDir(const char* dir) { fInputDir = (dir && dir[0] != '\0')
320 * Parse a string as a load option
322 * @param what String to pass
324 * @return Load option value, or 0 in case of error
326 static ETrees ParseLoad(const char* what);
329 @param o Object to copy from */
330 AliFMDInput(const AliFMDInput& o)
366 /** Assignement operator
367 @return REference to this */
368 AliFMDInput& operator=(const AliFMDInput&) { return *this; }
370 * Get user supplued data
379 virtual Float_t GetSignal(UShort_t d, Char_t r, UShort_t s, UShort_t t);
381 static const char* TreeName(ETrees tree, bool shortest=false);
384 * Make a chain of specified data
386 * @param what What data to chain. Possible values are
387 * - ESD Event summary data (AliESD)
388 * - MC Simulation data (galice)
389 * @param datadir Data directory to scan
390 * @param recursive Whether to recurse into sub-directories
392 * @return Pointer to newly create chain, or null
394 static TChain* MakeChain(const char* what, const char* datadir,
395 bool recursive=false);
397 * Scan a directory (optionally recursive) for data files to add to
398 * the chain. Only ROOT files, and files which name contain the
399 * passed pattern are considered.
401 * @param dir Directory to scan
402 * @param chain Chain to add data to
403 * @param pattern Pattern that the file name must contain
404 * @param recursive Whether to scan recursively
406 static void ScanDirectory(TSystemDirectory* dir,
407 const TString& olddir,
409 const char* pattern, bool recursive);
411 TString fGAliceFile; // File name of gAlice file
412 AliRunLoader* fLoader; // Loader of FMD data
413 AliRun* fRun; // Run information
414 AliStack* fStack; // Stack of particles
415 AliLoader* fFMDLoader; // Loader of FMD data
416 AliRawReader* fReader; // Raw data reader
417 AliFMDRawReader* fFMDReader; // FMD raw reader
418 AliFMD* fFMD; // FMD object
419 AliESDFMD* fESD; // FMD ESD data
420 AliESDEvent* fESDEvent; // ESD Event object.
421 TTree* fTreeE; // Header tree
422 TTree* fTreeH; // Hits tree
423 TTree* fTreeTR; // Track Reference tree
424 TTree* fTreeD; // Digit tree
425 TTree* fTreeS; // SDigit tree
426 TTree* fTreeR; // RecPoint tree
427 TTree* fTreeA; // Raw data tree
428 TChain* fChainE; // Chain of ESD's
429 TClonesArray* fArrayE; // Event info array
430 TClonesArray* fArrayH; // Hit info array
431 TClonesArray* fArrayTR; // Hit info array
432 TClonesArray* fArrayD; // Digit info array
433 TClonesArray* fArrayS; // SDigit info array
434 TClonesArray* fArrayR; // Rec points info array
435 TClonesArray* fArrayA; // Raw data (digits) info array
436 AliHeader* fHeader; // Header
437 TGeoManager* fGeoManager; // Geometry manager
438 Int_t fTreeMask; // Which tree's to load
439 TString fRawFile; // Raw input file
440 TString fInputDir; // Input directory
441 Bool_t fIsInit; // Have we been initialized
442 Int_t fEventCount; // Event counter
443 Int_t fNEvents; // The maximum number of events
444 static const ETrees fgkAllLoads[kUser+1]; // List of all possible loads
445 ClassDef(AliFMDInput,0) //Hits for detector FMD
448 inline Bool_t AliFMDInput::ProcessHit(AliFMDHit*,TParticle*) { return kTRUE; }
449 inline Bool_t AliFMDInput::ProcessTrackRef(AliTrackReference*,
450 TParticle*) { return kTRUE; }
451 inline Bool_t AliFMDInput::ProcessTrack(Int_t,TParticle*,
452 AliFMDHit*) { return kTRUE; }
453 inline Bool_t AliFMDInput::ProcessParticle(Int_t,TParticle*) { return kTRUE; }
454 inline Bool_t AliFMDInput::ProcessDigit(AliFMDDigit*) { return kTRUE; }
455 inline Bool_t AliFMDInput::ProcessSDigit(AliFMDSDigit*) { return kTRUE; }
456 inline Bool_t AliFMDInput::ProcessRawDigit(AliFMDDigit*) { return kTRUE; }
457 inline Bool_t AliFMDInput::ProcessRawCalibDigit(AliFMDDigit*) { return kTRUE; }
458 inline Bool_t AliFMDInput::ProcessRecPoint(AliFMDRecPoint*) { return kTRUE; }
459 inline Bool_t AliFMDInput::ProcessESD(UShort_t,Char_t,UShort_t,UShort_t,
460 Float_t,Float_t) { return kTRUE; }
461 inline Bool_t AliFMDInput::ProcessUser(UShort_t,Char_t,UShort_t,UShort_t,
462 Float_t) { return kTRUE; }
463 inline Float_t AliFMDInput::GetSignal(UShort_t, Char_t, UShort_t, UShort_t) {
468 //____________________________________________________________________